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FINAL  REPORT 

AFOSR  CONTRACT  NO  F49620-82-K-0033 

9/1/82  through  12/31/87 

I.  Experimental  Measurements  and  “ Bulk  Flow ”  Model  Developments  -  D.  Childs 

The  work  conducted  in  this  portion  of  the  research  contract  has  consisted  of  the 
following  tasks: 

(a)  the  development  of  a  facility  and  apparatus  for  measuring  the  leakage,  axial  pres¬ 
sure  profiles,  and  rotordynamic  (stiffness  and  damping)  coefficients  of  labyrinth 
seals, 

(b)  the  measurement  of  the  test  data  cited  above  for  a  range  of  labyrinth-seal  con¬ 
figurations,  and 

(c)  the  development  and  validation  of  “bulk-flow”  models  for  the  prediction  of  leakage 
and  rotordynamic  coefficients  of  labyrinth  seals. 

All  of  these  objectives  have  been  met  in  full  and  are  documented  in  the  following  journal 
publications: 

•  Childs,  D.  W.  and  Scharrer,  J.  K.,  “Experimental  Rotordynamic  Coefficient  and 
Results  for  Teeth-On-Rotor  and  Teeth-On-Stator  Labyrinth  Gas  Seals,”  ASME 
Trans.  J.  of  Engineering  for  Gas  Turbine  and  Power,  Vol.  108,  October  1986, 
pp.  599-604. 

•  Childs,  D.  W.,  Nelson,  C.  E.,  Nicks,  C.,  Scharrer,  J.,  Elrod,  D.,  and  Hale,  K., 
“Theory  Versus  Experiment  for  the  Rotordynamic  Coefficients  of  Annular  Gas 
Seals:  Part  1-Test  Facility  and  Apparatus”,  ASME  Transaction  Journal  of  Tri¬ 
bology,  July  1986,  Vol.  108,  pp.  426-432. 

•  Nelson,  C.,  Childs,  D.,  Nicks,  C.,  Elrod,  D.,  and  Hale,  K.,  “Theory  Versus  Exper¬ 
iment  for  the  Rotordynamic  Coefficients  of  Annular  Gas  Seals:  Part  2-Constant- 
Clearance  and  Convergent-  Tapered  Geometry,”  ASME  Transaction  Journal  of 
Tribology,  July  1986,  Vol.  108,  pp.  433-438. 

•  Childs,  D.,  and  Scharrer,  J.,  “An  Iwatsubo-Based  Solution  for  Labyrinth  Seals  - 
Comparison  to  Experimental  Results,”  ASME  Transaction  Journal  of  Engineer¬ 
ing  for  Gas  Turbine  and  Power,  April  1986,  Vol.  108,  pp.  325-331. 

•  Childs,  D.  and  Elrod,  D.,  “Rotordynamic  Coefficient  and  Leakage  Test  Results 
for  Interlock  and  Tooth-on- Stator  Labyrinth  Seals,”  accepted  for  presentation  at 
the  ASME  International  Gas  Turbine  Conference,  Amsterdam,  The  Netherlands, 
under  review  consideration  ASME  Trans,  for  Power. 


II.  Finite  Difference  Developments  and  Results  -  D.  Rhode 

This  portion  of  the  research  program  involved  the  following  tasks: 

(a)  the  development  of  a  3-D  finite  difference  computer  code,  using  a  generalized 
body-fitted  coordinate  system,  for  predicting  forces  as  well  as  the  distribution  of 
various  other  fiowfield  quantities 

(b)  the  development  of  a  procedure  for  calculating  the  rotordynamic  force  for  a  seal 
with  any  arbitrary  number  of  cavities  at  an  affordable  CPU  expenditure 

(c)  the  computation  of  forces  at  various  operating  conditions 

These  objectives  have  been  met  and  are  documented  in  the  following: 

•  Rhode,  D.  L.  and  Hensel,  S.  J.,  “Three-Dimensional  Computation  of  Rotordy¬ 
namic  Force  Distributions  in  a  Labyrinth  Seal,”  accepted  for  presentation  at  the 
ASME/AIAA  First  National  Fluid  Dynamics  Congress,  Cincinnati,  OH,  24-28 
July  1988. 

•  Rhode,  D.  L.  and  Hensel,  S.  J.,  “Labyrinth  Seal  Rotordynamic  Forces  Predicted 
with  a  Three-Dimensional  Navier-Stokes  Computer  Code,”  accepted  for  presentar 
tion  at  the  24th  ALAA/ASEE/ASME/SAE  Joint  Propulsion  Conference,  Boston, 
MA,  11-14  July  1988. 

•  Rhode,  D.  L.  and  Nail,  G.  H.,  “Computation  of  Cavity-By-Cavity  Flow  Devel¬ 
opment  in  Generic  Labyrinth  Seals,”  submitted  for  presentation  at  the  ASME 
International  Computers  in  Engineering  Conference,  San  Francisco,  CA,  31  July 

-  3  August  1988. 

•  Rhode,  D.  L.  and  Sobolik,  S.  R.,  “Simulation  of  Subsonic  Flow  Through  a  Generic 
Labyrinth  Seal,”  ASME  Trans.  Journal  of  Engineering  for  Gas  Turbines  and 
Power,  October  1986,  Vol.  108,  pp.  674-680. 

•  Rhode,  D.  L.,  Dernko,  J.  A.,  Traegner,  U.  K.,  Morrison,  G.  L.  and  Sobolik,  S.  R., 
“Prediction  of  Incompressible  Flow  in  Labyrinth  Seals,”  ASME  Trans.  Journal 
of  Fluids  Engineering,  March  1986,  Vol.  108,  pp.  19-25. 

Appendix  B  contains  three  of  the  most  recent  papers: 

•  Rhode,  D.  L.  and  Hensel,  S.  J.,  “Three-Dimensional  Computation  of  Rotordy¬ 
namic  Force  Distributions  in  a  Labyrinth  Seal,”  accepted  for  presentation  at  the 
ASME/AIAA  First  National  Fluid  Dynamics  Congress,  Cincinnati,  OH,  24-28 
July  1988. 

•  Rhode,  D.  L.  and  Hensel,  S.  J.,  “Labyrinth  Seal  Rotordynamic  Forces  Predicted 
with  a  Three-Dimensional  Navier-Stokes  Computer  Code,”  accepted  for  presenta¬ 
tion  at  the  24th  AIAA/ASEE/ASME/SAE  Joint  Propulsion  Conference,  Boston, 
11-14  July  1988. 

•  Rhode,  D.  L.  and  Nail,  G.  H.,  “Computation  of  Cavity-By-Cavity  Flow  Devel¬ 
opment  in  Generic  Labyrinth  Seals,”  submitted  for  presentation  at  the  ASME 
International  Computers  in  Engineering  Conference,  San  Francisco,  CA,  31  July 

-  3  August  1988. 


The  following  students  have  participated  in  this  portion  of  the  AFRAPT  program: 
Steve  Sobolik,  Steve  Hensel,  Greg  Nail,  and  Robert  Hibbs.  Steve  Hensel  worked  at  Garrett, 
one  summer,  finished  his  M.  S.  Thesis  in  December  1987  and  is  continuing  his  studies  for  a 
Ph.D.  A  summer  position  was  not  available  in  time  for  Steve  Sobolik  and  Greg  Nail.  Steve 
finished  his  M.  S.  Thesis  in  August  1984  and  went  to  work  for  Sandia  National  Labs.  Greg 
finished  his  M.  S.  Thesis  in  December  1987  and  is  continuing  his  education  for  the  Ph.D. 
degree.  Robert  Hibbs  worked  last  summer  at  United  Technologies  in  East  Hartford  and 
will  soon  finish  his  M.  S.  Thesis. 

We  are  grateful  for  the  sponsorship  of  this  research  program  by  AFOSR.  As  with  the 
test  facility,  the  computational  capability  continues  to  be  unique  throughout  the  world. 
The  value  of  this  widely  applicable  predictive  capability  has  clearly  been  demonstrated. 
Also,  previously  unavailable  values  of  various  shear  stress  quantities  have  been  computed, 
which  were  particularly  useful  in  refining  the  bulk-flow  model  mentioned  in  the  previous 
section. 

Further  work  is  important  and  should  include  an  analysis  of  the  effect  of  geometric 
configuration  including  stepped  seals,  as  well  as  the  effect  of  rotor  orbits  whose  centers  are 
eccentric  with  respect  to  the  housing. 


Copies  of  the  three  most-recent  publications  are  included  in  Appendix  A. 

•  Hawkins,  L.,  Childs,  D.,  and  Hale,  K.,  “Experimental  Results  for  Labyrinth  Gas 
Seals  with  Honeycomb  Stators:  Comparisons  to  Smooth-Stator  Seals  and  Theo¬ 
retical  Predictions,”  submitted  for  the  1988  ASME-ASLE  Joint  Tribology  Con¬ 
ference  and  ASME  Journal  of  Tribology. 

•  Scharrer,  J.,  “Theory  versus  Experiment  for  the  Rotordynamic  Coefficients  of 
Labyrinth  Gas  Seals:  Part  I  -  A  Two  Control  Volume  Model,”  Rotating  Machinery 
Dynamics,  Vol.  2,  ASME  1987,  pp.  427-434,  accepted  for  presentation,  ASME 
Journal  of  Vibration,  Acoustics,  Stress,  and  Reliability  in  Design. 

•  Childs,  D.  and  Scharrer,  J.,  “Theory  Versus  Experiment  for  the  Rotordynamic 
Coefficients  of  Labyrinth  Gas  Seals:  Part  H  -  A  Comparison  to  Experiment,” 
Rotating  Machinery  Dynamics,  Vol.  2,  ASME  1987,  pp.  427-434,  accepted  for 
presentation,  ASME  Journal  of  Vibration,  Acoustics,  Stress,  and  Reliability  in 
Design. 

The  AFRAPT  participation  in  this  program  has  included  the  students:  Joseph  Schar¬ 
rer  and  Lawrence  Hawkins.  Joe  worked  at  G.E.  Lynn  during  the  summer,  completed  his 
Ph.D.  in  January  1987  and  is  continuing  to  work  in  rotating  machinery  with  Rocketdyne. 
Larry  worked  for  Garrett  during  two  summers,  completed  his  M.S.  degree  in  January  1988, 
and  is  also  going  to  work  for  Rocketdyne. 

The  support  provided  by  AFOSR  is  deeply  appreciated.  The  test  apparatus  we  have 
developed  with  AFOSR  support  continues  to  be  unique  in  the  world.  The  test  results 
and  models  have  been  of  extraordinary  value  in  resolving  rotordynamics  stability  issues 
and  have  been  used  directly  to  develop  higher-performance  commercial  turbomachinery. 
Considerable  work  remains  to  be  done!  No  data  or  models  are  available  yet  for  stepped 
seals,  see-through  labyrinth  seals  at  reduced  L/D  ratios,  brush  seals,  etc. 


Appendix  A 


Hawkins,  L.,  Childs,  D.,  and  Hale,  K.,  “Experimental  Results  for  Labyrinth  Gas  Seals  with 
Honeycomb  Stators:  Comparisons  to  Smooth-Stator  Seals  and  Theoretical  Predictions," 
submitted  for  the  1988  ASME-ASLE  Joint  Tribology  Conference  and  ASME  Journal  of 
Tribology. 

Scharrer,  J.,  “Theory  versus  Experiment  for  the  Rotordynamic  Coefficients  of  Labyrinth 
Gas  Seals:  Part  I  -  A  Two  Control  Volume  Model,”  Rotating  Machinery  Dynamics,  Vol.  2, 
ASME  1987,  pp.  427-434,  accepted  for  presentation,  ASME  Journal  of  Vibration,  Acous¬ 
tics,  Stress,  and  Reliability  in  Design. 

Childs,  D.  and  Scharrer,  J.,  “Theory  Versus  Experiment  for  the  Rotordynamic  Coefficients 
of  Labyrinth  Gas  Seals:  Part  H  -  A  Comparison  to  Experiment,”  Rotating  Machinery 
Dynamics,  Vol.  2,  ASME  1987,  pp.  427-434,  accepted  for  presentation,  ASME  Journal  of 
Vibration,  Acoustics,  Stress,  and  Reliability  in  Design. 
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EXPERIMENTAL  RESULTS  FOR  LABYRINTH  GAS  SEALS 
WITH  HONEYCOMB  STATORS:  COMPARISONS  TO 
SMOOTH-STATOR  SEALS  AND  THEORETICAL  PREDICTIONS1 

Larry  Hawkins,  Member  Technical  Staff 
Rocketdyne  Division,  Rockwell  International 
Canoga  Park,  CA  91304 

Dara  Childs,  Professor 
Keith  Hale,  Research  Engineer 
Texas  A&M  University 
College  Station,  Texas  77843 


ABSTRACT 

Experimental  measurements  are  presented  for  the  rotordynamic  stiffness  and  damping 
coefficients  of  a  teeth-on-rotor  labyrinth  seal  with  a  honeycomb  stator .  Inlet  circumferential 
velocity,  inlet  pressure,  rotor  speed,  and  seal  clearance  are  primary  variables.  Results  are 
compared  to  (a)  data  for  teeth-on-rotor  labyrinth  seals  with  smooth  stators,  and  (b)  ana¬ 
lytical  predictions  from  a  two-control-volume  compressible  flow  model.  The  experimental 
results  show  that  the  honeycomb-stator  configuration  is  more  stable  than  the  smooth-stator 
configuration  at  low  rotor  speeds.  At  high  rotor  speeds,  the  stator  surface  does  not  affect 
stability.  The  theoretical  model  predicts  the  cross-coupled  stiffness  of  the  honeycomb-stator 
seal  correctly  within  25%  of  measured  values.  The  model  provides  accurate  predictions  of 
direct  damping  for  large  clearance  seals;  however,  the  model  predictions  and  test  results 
diverge  with  increasing  running  speed.  Overall,  the  model  does  not  perform  as  well  for  low 
clearance  seals  as  for  high  clearance  seals. 

INTRODUCTION 

Modern  turbomachines  can  be  subject  to  problems  due  to  unstable  or  self-excited 
motion.  This  type  of  motion  is  characterized  by  a  rotor  whirling  at  a  natural  frequency  ' 
that  is  less  than  ito  rotational  speed.  The  unstable  motion  is  caused  by  a  tangential 
force  acting  on  the  rotor  in  its  whirl  direction.  Several  excitation  mechanisms  have  been 
identified  for  unstable  motion  (Ehrich  and  Childs,  1984);  one  of  these  is  the  force  developed 
in  a  labyrinth  seal. 

^his  work  was  supported  in  part  by  NASA  Grant  NAG3-181  from  NASA  Lewis  Research 
Center  (Technical  Monitor,  Robert  Hendricks)  and  AFOSR  Contract  F49620-82-K-0033 
(Technical  Monitor,  Tony  Amos) 
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Labyrinth  seal  forces  are  characterized  by  rotordynamic  stiffness  and  damping  coef¬ 
ficients.  For  small  motion  of  a  seal  rotor  about  a  centered  position,  the  rotordynamic 
coefficients  are  defined  by  the  following  force-motion  model 


(1) 


where  X  and  Y  define  the  motion  of  the  seal  rotor  relative  to  the  seal  stator,  and  Fx 
and  Fy  are  the  reaction  force  components  acting  on  the  seal  rotor.  The  rotordynamic 
coefficients  ( K ,  k ,  C,  and  c)  represent  the  direct  stiffness,  cross-coupled  stiffness,  direct 
damping,  and  cross-coupled  damping  respectively. 

This  report  presents  experimental  measurements  of  the  rotordynamic  coefficients  for  a 
teeth-on-rotor  labyrinth  seal  with  a  honeycomb  stator.  The  first  systematic  test  program 
for  measuring  seal  rotordynamic  coefficients  was  performed  at  the  Technical  University 
of  Stuttgart  (Benckert  and  Wachter, 1978, 1980), (Benckert, 1980).  Stiffness  data  were  pub¬ 
lished  for  three  types  of  seals:  (a)  teeth-on-stator,  (b)  teeth-on-rotor  and  stator,  and 
(c)  teeth-on-stator  and  steps  or  grooves  on  the  rotor.  Wright  (1983)  published  data  on 
equivalent  radial  and  tangential  stiffnesses  for  single-cavity,  teeth-on-stator  seals.  Childs 
and  Scharrer  (1986,1987)  investigated  teeth-on-rotor  and  teeth-on-stator  labyrinth  seals  at 
Texas  A&M  University.  They  measured  stiffness  and  damping  coefficients  while  varying  in¬ 
let  tangential  velocity,  rotor  speed,  inlet  pressure,  and  clearance.  Elrod  and  Childs  (1988) 
investigated  smooth-rotor/honeycomb-stator  annular  seals,  and  compared  the  results  to 
data  for  smooth-rotor/smooth-stator  annular  seals.  They  found  the  honeycomb-stator 
seals  to  be  more  stable  than  smooth-stator  seals  throughout  the  range  of  variables  tested. 

The  labyrinth-rotor/honeycomb-stator  configuration  was  tested  for  several  reasons: 
(a)  This  combination  is  a  common  industrial  configuration,  (b)  no  test  data  for  this  com¬ 
bination  exists  in  the  published  literature,  and  (c)  the  results  of  Elrod  and  Childs  (1988) 
indicate  that  seals  with  honeycomb  stators  may  have  a  stability  advantage  over  smooth- 
stator  seals.  The  last  item  motivates  the  comparison  in  this  report  of  the  honeycomb-stator 
seal  data  to  the  data  of  Childs  and  Scharrer  (1987)  for  labyrinth-rotor/smooth-stator  seals. 
The  honeycomb-stator  seal  data  are  also  compared  to  the  predictions  of  Scharrer’s  (1987) 
theoretical  model. 


THEORY 
Stability  Analysts 

As  a  consequence  of  equation  (1),  the  radial  and  tangential  forces  acting  on  a  seal 
rotor,  which  is  executing  a  circular  orbit  of  amplitude  A,  are  defined  by 

Fr  =  —{K  +  cu)A 
Ft  =  (k  -  Coj)A 


where  u>  is  the  running  speed.  Thus,  K  and  c  produce  a  radial  force  on  the  seal  rotor,  and 
k  and  C  produce,  a  tangential  force.  The  radial  force  Fr  in  a  labyrinth  seal  is  normally 
small  compared  to  other  radial  fores  acting  on  a  rotor;  therefore  ,  K  and  c  are  of  minor 
consequence.  If  positive,  the  tangential  force  Ft  is  destabilizing  since  it  supports  the 
whirling  motion  of  a  forward  whirling  rotor.  Both  k  and  C  are  positive  for  most  practical 
labyrinth  seal  applications;  hence,  the  compelling  reason  for  determining  the  rotordynamic 
seal  coefficients  is  to  determine  the  relative  values  of  k  and  C.  The  whirl  frequency  ratio, 
defined  by 

/  =  k/Cu, 

is  used  in  this  report  to  compare  k  and  C .  Reducing  the  value  of  /  improves  the  stability 
of  a  rotor  system. 

Scharrer’s  Analysis 

Prediction  of  the  rotordynamic  coefficients  requires  a  flow  field  model  to  predict  the  ax¬ 
ial  and  circumferential  pressure  distribution  acting  on  the  seal  rotor.  Most  early  attempts 
to  model  the  flow  field  in  a  labyrinth  seal  used  a  single  control  volume,  concentrating  on 
the  circumferential  flow  components.  However,  the  labyrinth  seal  is  known  to  have  two 
distinct  flow  regimes:  a  jet  flow  region  in  the  leakage  path  and  a  recirculation  region  in  the 
cavity.  Hence,  Fujikawa  et  al.  (1984),  Wyssmann  et  al.  (1984),  and  Scharrer  (1987)  have 
developed  two-control-volume  analyses  to  accurately  model  the  known  physics  of  the  flow. 
Scharrer’s  model  is  used  in  this  report  to  generate  theoretical  predictions  for  comparison 
to  experimental  data. 

In  Scharrer’s  model,  the  shear  stresses  at  the  rotor  and  stator  surfaces  (r#  and  rr)  are 
modeled  using  the  Blasius  formula  for  turbulent  pipe  flow 


where  Um  is  the  mean  flow  velocity  relative  to  the  surface  upon  which  the  shear  stress  is 
cting,  and  Dh  is  the  hydraulic  diameter  of  the  particular  control  volume.  No  published 
data  are  available  for  the  friction  coefficients,  no  and  mo,  for  the  honeycomb-stator  surface 
used  in  the  *ests  reported  here.  The  following  values  were  determined  empirically  from 
pressure  drop  versus  flow  tests  for  smooth-rotor/honeycomb-stator  seals 

ms  =  —0.1083  ns  =  0.2820. 

Scharrer  uses  a  perturbation  analysis  to  linearize  the  governing  equations.  Expanding 
the  governing  equations  in  the  perturbation  variables  yields  a  system  of  twelve  linear 
algebraic  equations  per  cavity.  Solution  of  these  equations  yields  the  pressure  distribution 
along  and  around  the  seal.  Integration  of  the  pressure  distribution  leads  to  the  solution 
for  the  rotordynamic  coefficients. 
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TEST  APPARATUS 

The  test  results  reported  here  were  obtained  using  the  Texas  A&M  Air  Seal  Test 
Apparatus.  The  test  apparatus  will  be  described  here  briefly.  A  detailed  description  was 
provided  by  Childs  et  al.  (1986).  As  illustrated  in  figure  1,  the  rotor  shaft  is  suspended, 
pendulum  fashion,  from  an  upper,  rigidly-mounted,  pivot  shaft.  This  arrangement  allows 
horizontal  (harmonic)  motion  of  the  rotor.  A  cam  within  the  pivot  shaft  provides  vertical 
(static)  positioning  of  the  rotor.  The  rotor  is  excited,  horizontally,  by  a  hydraulic-shaker 
head  which  acts  on  the  rotor-shaft  housing.  The  design  of  the  test  rig,  which  is  further 
illustrated  in  figure  2,  permits  the  installation  of  various  rotor/stator  combinations.  The 
test  apparatus  has  been  modified  since  the  1986  reference  to- permit  an  increase  in  top 
operating  speeds  from  8,000  to  16,000  cpm.  Changes  include  the  use  of  a  hydraulically 
fitted  rotor,  the  introduction  of  high-speed  carbon  seals,  and  the  replacement  of  a  roller- 
element  thrust  bearing  with  a  Torrington,  water-iubricated,  swing-pad  bearing.  The  stator 
of  figure  2  is  supported  in  the  test  section  housing  by  three  piezo-electric,  quartz,  load  cells 
in  a  trihedral  configuration.  These  load  cells  measure  the  pressure-induced  forces  due  to 
rotor  motion  within  the  stator.  Accelerometers  are  provided  on  the  stator  to  correct  for 
acceleration-induced  forces  which  are  measured  by  the  load  cell. 

TEST  PARAMETERS 

The  parameters  varied  during  the  tests  were  supply  pressure,  rotor  speed,  circumfer¬ 
ential  velocity  of  the  inlet  air,  and  seal  configuration.  Test  results  for  six  teeth-on-rotor 
labyrinth  gas  seal  configurations  are  presented.  Three  of  the  seals  have  honeycomb  stators, 
each  with  a  different  rotor-to-stator  clearance.  The  other  three  seals  have  smooth  stators, 
each  with  a  clearance  equal  to  one  of  the  honeycomb-stator  seals.  The  seals  are  illustrated 
in  figure  3.  Seals  1,  2,  3,  and  4  were  tested  for  this  study,  and  the  data  for  these  seals  are 
reported  here  for  the  first  time.  Seals  5  and  6  were  tested  previously  and  documented  by 
Childs  and  Scharrer  (1987).  The  data  are  presented  here  again  to  provide  comparison  to 
the  corresponding  honeycomb  stator  seals  (seals  2  and  3). 

The  test  points  for  the  other  three  variables  are  shown  in  table  1.  The  inlet  air 
pressure  and  attendant  mass  flow  rate  through  the  seal  are  controlled  with  a  flow  control 
valve  located  upstream  of  the  seal.  Rotor  speed  is  controlled  by  a  variable  speed  electric 
motor  with  a  frequency  controller.  The  inlet  circumferential  velocities  are  controlled  using 
the  inlet  guide  vanes  shown  in  figure  4.  The  guide  vanes  are  contained  in  sleeves  and  located 
immediately  upstream  of  the  test  seal.  The  no-prerotation  case  is  obtained  without  guide 
vanes.  High  and  low  prerotation  velocities  are  obtained  for  the  different,  guide-vane-depths 
A  of  figure  4.  The  inlet  circumferential  velocity  is  calculated  from  measured  values  for  the 
volumetric  flow  rate,  upstream  temperature  and  pressure,  and  a  flow-turning  correction 
in  accordance  with  Cohen  et  al.  (1972).  The  circumferential  velocity  can  not  be  varied 
arbitrarily,  because  it  depends  on  the  supply  pressure  and  the  flow  resistance  of  the  seal 
being  tested.  The  inlet  circumferential  velocity  set  points  for  each  seal  are  shown  in  table 
2.  Inlet  circumferential  velocity  varies  a  maximum  of  ±  10%  with  supply  pressure  and 
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rotor  speed  from  the  values  in  table  2.  When  reviewing  the  following  test  results,  figure 
3,  table  1,  and  table  2  should  be  consulted  for  the  definitions  of  symbols  used.  Also,  note 
that  circumferential  velocity  ratio  is  used  to  represent  inlet  circumferential  velocity;  this 
is  the  ratio  of  inlet  circumferential  velocity  to  rotor  surface  speed. 

As  noted  in  table  2,  the  inlet  circumferential  velocity  test  point  varies  significantly 
when  seal  clearance  is  varied.  Therefore,  the  rotordynamic  coefficients  cannot  be  plotted 
versus  clearance  because  the  variation  in  inlet  circumferential  velocity  with  clearance  would 
afTect  the  results.  The  effect  of  clearance  is  displayed  by  plotting  the  coefficients  versus 
inlet  circumferential  velocity  for  each  seal  on  the  same  plot.  This  procedure  allows  only 
one  rotor  speed  and  one  supply  pressure  per  plot. 

The  uncertainty  of  the  rotordynamic  coefficients  was  calculated  using  the  method 
described  by  Holman  (1978).  For  the  seal  configurations  tested,  the  maximum  uncertainties 
in  the  stiffness  and  damping  coefficients  were  15  N/mm  (86  Ib/in)  and  32  N-s/m  (0.18  1b- 
s/in),  respectively.  The  uncertainty  in  the  cross-coupled  damping  coefficients  are  of  the 
same  order  of  magnitude  as  the  coefficients  themselves;  therefore,  cross-coupled  damping 
data  are  not  presented  here. 

TEST  RESULTS 

Honeycomb-stator  seal  data  are  compared  to  smooth-stator  seal  data  in  figures  5-12. 
In  these  figures,  solid  lines  represent  the  honeycomb-stator  test  results  and  broken  lines 
represent  the  smooth-stator  test  results.  The  smallest  clearance  seals  (seals  1  &  4)  and  the 
highest  inlet  circumferential  velocity  (swirl  3)  are  used  for  plots  in  which  these  parameters 
.are  not  varied.  Further  data  are  presented  by  Hawkins  (1988). 

Leakage 

Leakage  is  represented  by  the  flow  coefficient, 

mjm 

nDCrP, ' 

Figure  5  is  a  plot  of  flow  coefficient  versus  seal  clearance  for  an  inlet  pressure  of  3.08  bar 
and  a  rotor  speed  of  3000  cpm.  Curve  1  represents  the  honeycomb-stator  seals  and  curve  4 
represents  the  smooth  stator  seals.  Leakage  did  not  vary  with  inlet  circumferential  velocity, 
thus  the  data  presented  are  for  inlet  circumferential  velocity  3  only.  Examination  of  the 
figure  reveals  that  the  honeycomb-stator  seal  leaks  more  at  the  smallest  clearance  and 
the  smooth-stator  seal  leaks  more  at  the  largest  clearance.  This  result  is  consistent  with 
Stocker  et  al.  (1977).  One  would  expect  the  honeycomb-stator  seal  to  leak  less  than  the 
smooth-stator  seal  because  its  greater  roughness  increases  the  flow  resistance.  However, 
in  the  honeycomb-stator  seal,  part  of  the  leakage  air  may  bypass  the  apparent  flow  path 
by  passing  into  and  out  of  honeycomb  cells.  Thus,  relative  to  the  smooth-stator  seal,  the 
honeycomb-stator  6eal  has  a  larger  effective  leakage  area  for  a  given  clearance.  When  the 
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clearance  is  small,  the  increased  leakage  area  results  in  greater  leakage  in  the  honeycomb- 
stator  seal  versus. the  smooth-stator  seal. 

Figure  6  is  similar  to  figure  5,  except  that  leakage  is  represented  by  the  dimensional 
mass  flow  rate.  This  figure  shows  that  leakage  increases  as  clearance  increases. 

Direct  Stiffness 

Direct  stiffness  is  plotted  versus  rotor  speed  for  various  pressures  in  figure  7.  Again, 
the  solid  lines  are  for  the  smallest  clearance  honeycomb-stator  seal,  and  the  broken  lines 
are  for  the  smallest  clearance  smooth-stator  seal.  Direct  stiffness  of  the  honeycomb-stator 
seal  is  negative  and  becomes  more  negative  as  rotor  speed  increases.  The  smooth-stator 
seal  has  a  similar  characteristic,  but  has  a  larger  direct  stiffness  magnitude.  Direct  stiffness 
becomes  more  negative  as  pressure  increases  for  both  stator  surfaces. 

Cross-Coupled  Stiffness 

Cross-coupled  stiffness  is  plotted  versus  rotor  speed  for  various  pressures  in  figure  8. 
Cross-coupled  stiffness  increases  with  rotor  speed  for  both  seals.  For  the  honeycomb-stator 
seal,  cross-coupled  stiffness  is  negative  at  low  speed  and  increases  to  about  300  N/mm  at 
the  highest  rotor  speed.  For  the  smooth  stator  seal,  cross-coupled  stiffness  has  a  small 
positive  value  at  low  rotor  speeds,  increasing  to  about  350  N/mm  at  the  highest  rotor 
speed.  Due  to  the  results  of  Elrod  and  Childs  (1988),  cross-coupled  stiffness  was  expected 
to  be  less  positive  for  the  honeycomb-stator  seal  compared  to  the  smooth-stator  seal  for  all 
rotor  speeds.  The  data  show  that  the  cross-coupled  stiffness  of  the  two  seals  have  similar 
magnitudes  at  high  rotor  speeds.  This  plot  also  shows  that  increasing  pressure  increases 
the  absolute  value  of  cross-coupled  stiffness. 

Figure  9  illustrates  cross-coupled  stiffness  versus  circumferential  velocity  for  the  three 
honeycomb-stator  seals  (1,  2,  3)  and  the  three  smooth-stator  seals  (4,  5,  6).  Figure  9(a) 
is  for  an  inlet  pressure  of  3.08  bar  and  rotor  speed  of  3000  cpm,  and  figure  9(b)  is  for 
the  same  pressure  and  rotor  speed  of  16,000  cpm.  Figure  9(a)  shows  that  cross-coupled 
stiffness  increases  roughly  linearly  as  inlet  circumferential  ratio  increases  for  both  stator 
surfaces.  The  increase  is  greater  with  the  smooth-stator  seal.  In  figure  9(b),  cross-coupled 
stiffness  increases  significantly  from  zero  inlet  circumferential  velocity  to  the  first  positive 
value  of  inlet  velocity;  however,  from  the  first  positive  inlet  velocity  to  the  second  positive 
inlet  velocity,  the  cross-coupled  stiffness  increases  only  slightly  or  in  some  cases  decreases. 
This  trend  is  present  for  both  stator  surfaces.  Returning  to  figure  9(a),  for  the  honeycomb- 
stator  seal,  cross-coupled  stiffness  is  much  lower  in  the  small  clearance  seal  (seal  1)  than 
in  the  larger  clearance  seals  (seals  2  Ac  3).  Increasing  seal  clearance  has  little  effect  on  the 
cross-coupled  stiffness  of  the  smooth-stator  seal.  At  the  higher  rotor  speed  (figure  9(b)), 
seal  clearance  has  little  effect  on  cross-coupled  stiffness  in  the  honeycomb-stator  seal.  The 
small  clearance  smooth-stator  seal  (seal  4)  has  a  much  higher  cross-coupled  stiffness  than 
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the  larger  clearance  seals  (seals  5  <Sc  6).  Note  that  cross-coupled  stiffness  of  the  larger 
clearance  smoothrstator  seals  does  not  depend  on  rotor  speed. 

Direct  Damping 

Direct  damping  is  plotted  versus  rotor  speed  for  various  pressures  in  figure  10.  Direct 
damping  has  essentially  the  same  magnitude  for  either  stator  surface.  However,  damp¬ 
ing  for  the  honeycomb-stator  seal  first  increases  and  then  decreases  with  increasing  rotor 
speed,  while  damping  in  the  smooth-stator  seal  does  not  depend  on  rotor  speed.  Also, 
damping  increases  with  pressure  for  both  stator  surfaces,  but  the  increase  is  greater  in  the 
honeycomb-stator  seal. 

Figure  11  is  a  plot  of  direct  damping  versus  inlet  circumferential  velocity  ratio  for  each 
seal  configuration.  Damping  does  not  have  a  clear  dependence  on  inlet  circumferential 
velocity  in  the  honeycomb-stator  seals,  but  increases  with  increasing  inlet  circumferential 
velocity  in  the  smooth-stator  seals.  Damping  increases  somewhat  from  seal  1  (the  smallest 
clearance  seal)  to  seal  2.  However,  damping  in  seal  2  and  seal  3  is  roughly  the  same.  In 
the  smooth-stator  seal,  damping  first  falls  and  then  rises  as  clearance  increases.  This  plot 
is  for  pressure  of  3.08  bar  and  rotor  speed  of  3000  cpm.  The  same  trends  are  repeated  at 
other  inlet  pressures  and  rotor  speeds. 

Whirl  Frequency  Ratio 

Figure  12  is  a  plot  of  whirl  frequency  ratio  versus  rotor  speed  for  the  three  honeycomb- 
stator  seals  and  the  three  smooth-stator  seals.  As  noted  previously,  lower  values  of  whirl 
frequency  ratio  indicate  a  more  stable  seal.  Note  that  at  the  lower  rotor  speeds,  each 
honeycomb-stator  seal  is  more  stable  than  the  comparable  smooth-stator  seal,  and  that 
for  a  particular  stator  surface,  the  small  clearance  seals  are  more  stable  than  the  larger 
clearance  seals.  At  the  higher  rotor  speeds,  stator  surface  does  not  significantly  affect 
stability.  Also,  the  small  clearance  seals  (seals  1  &  4)  are  less  stable  than  the  others  at  the 
higher  speeds. 

COMPARISON  TO  THEORY 

Data  for  the  honeycomb-stator  seals  are  compared  to  the  predictions  of  Scharrer’s 
model  in  figures  13-16.  The  experimental  data  are  represented  by  solid  lines  and  theoretical 
predictions  are  represented  by  broken  lines  in  each  figure. 

Cross-Coupled  Stiffness 
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Figure  13  is  a  plot  of  cross-coupled  stiffness  versus  rotor  speed  for  various  pressures 
for  the  smallest  clearance  seal.  The  theory  correctly  predicts  that  cross-coupled  stiffness 
rises  as  rotor  speed  rises.  The  theory  does  not  predict  the  negative  values  found  at  low 
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rotor  speeds.  In  general,  the  model  predicts  positive  values  of  cross-coupled  stiffness  cor¬ 
rectly  within  25%.  Figure  14  is  a  plot  of  cross-ccuplcd  stiffness  versus  inlet  circumferential 
velocity  ratio  for  the  three  honeycomb-stator  seals.  The  model  predicts  a  larger  rise  in 
cross-coupled  stiffness  with  increasing  circumferential  velocity  than  is  shown  in  the  ex¬ 
perimental  data.  The  model  also  predicts  little  increase  in  cross-coupled  stiffness  with 
increasing  clearance. 


Direct  Damping 

Direct  damping  is  plotted  versus  rotor  speed  for  various  pressures  in  figure  15.  The 
model  predicts  that  damping  increases  with  increasing  rotor  speed,  whereas  the  data  show 
that  damping  first  increases  and  then  decreases  with  increasing  rotor  speed.  The  model 
consistently  underpredicts  the  magnitude  of  the  damping  at  low  rotor  speeds,  and  is  accu¬ 
rate  at  high  rotor  speeds.  However,  at  the  highest  rotor  speeds  tested,  the  model  predic¬ 
tions  and  the  test  results  are  diverging.  Figure  16  is  a  plot  of  direct  damping  versus  inlet 
circumferential  velocity  ratio  for  the  three  honeycomb-stator  seals.  The  model  correctly 
predicts  that  damping  does  not  depend  on  inlet  circumferential  velocity  ratio.  The  model 
also  predicts  that  damping  does  not  increase  significantly  as  clearance  increases.  This  is 
true  only  at  the  larger  clearances. 

CONCLUSIONS 


The  test  data  support  the  following  conclusions  for  the  labyrinth-ro tor/honey comb- 
stator  seals: 

1)  Cross-coupled  stiffness  is  generally  positive.  Cross-coupled  stiffness  increases  with 
rotor  speed  and  with  inlet  circumferential  velocity.  At  the  lower  rotor  speeds,  cross- 
coupled  stiffness  is  much  lower  for  the  smallest  clearance  6eal  than  for  the  other  two 
seals.  At  the  higher  rotor  speeds,  cross-coupled  stiffness  is  approximately  the  same 
value  regardless  of  clearance. 

2)  Direct  damping  is  positive,  and  is  much  lower  in  the  smallest  clearance  seal  than  in 
the  two  larger  clearance  seals. 

By  comparing  the  results  for  the  honeycomb-stator  and  smooth-stator  seals,  the  fol¬ 
lowing  conclusions  may  be  drawn: 

1)  The  honeycomb-stator  seals  leak  more  than  the  smooth-stator  seals  when  the  clearance 
is  small.  The  honeycomb-stator  seals  leak  less  when  the  clearance  is  large. 

2)  The  honeycomb-stator  seal  is  more  stable  at  low  rotor  speeds.  For  high  rotor  speeds 
stator  surface  does  not  affect  stability. 
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By  comparison  of  experimental  results  and  theoretical  predictions  for  honeycomb- 
stator  seals,  the  following  conclusions  may  be  drawn: 

1)  The  model  does  not  predict  the  negative  values  measured  for  cross-coupled  stiffness 
at  low  rotor  speeds. 

2)  The  model  consistently  predicts  the  positive  values  of  cross-coupled  stiffness  of  the 
honeycomb-stator  seal  correctly  within  25%  of  the  measured  values.  The  model  cor¬ 
rectly  predicts  the  weak  dependence  of  cross-coupled  stiffness  on  clearance  for  the 
larger  clearances. 

3)  The  model  incorrectly  predicts  that  direct  damping  increases  with  speed,  and  does 
not  predict  the  decrease  in  damping  at  small  clearance.  For  the  two  larger  clearance 
seals  the  model  produces  good  results  for  tested  rotor  speeds  above  12,000  cpm.  Below 
12,000  cpm,  the  model  underpredicts  direct  damping  by  50%. 

In  general,  Scharrer’s  model  gives  useful  results  for  cross-coupled  stiffness  in  the 
labyrinth-rotor/honeycomb-stator  seal  for  the  range  of  variables  tested.  Scharrer’s  model 
can  give  useful  results  for  direct  damping  in  the  labyrinth-rotor/honeycomb-stator  seal  by 
applying  a  rotor-speed  correction  factor  to  the  predicted  damping.  Overall,  the  model 
produces  better  results  for  the  larger  clearances.  The  increased  significance  of  unmodeled 
effects  at  smaller  clearances  -  such  as  the  unmodeled  leakage  path  through  the  honey¬ 
comb  cells  -  is  probably  responsible  for  the  reduced  performance  of  the  model  at  smaller 
clearances. 

Values  of  the  rotordynamic  coefficients  for  the  two  larger  clearance  seals  tend  to  be 
much  closer  together  than  to  the  smaller  clearance  seal.  This  is  true  for  both  the  labyrinth- 
rotor/honeycomb-stator  seal  and  for  the  previously  untested  smallest  clearance  labyrinth- 
rotor/ smooth-stator  seal.  Since  there  are  many  practical  applications  where  labyrinth  seals 
are  used  with  clearances  below  the  tested  range,  further  testing  with  smaller  clearances 
are  required. 
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NOMENCLATURE 


'ms,  ns 
P 
R 
T 

X,Y 


Direct  and  cross-coupled  damping  coefficients  ( Ft/L ) 
Radial  clearance  (Z>) 

Rotor  diameter  (L) 

Seal  reaction-force  magnitude  ( F ) 

Direct  and  cross-coupled  stiffness  coefficients  (F / L) 
Whirl  frequency  ratio 
Leakage  mass  flow  rate  (M/Lt) 

Friction  coefficients 
Fluid  pressure  (F/L2) 

Gas  constant  for  air  (L2/Tf2) 

Fluid  temperature  ( T ) 

Rotor-to-stator  relative  displacement  components  ( L ) 
Kinematic  viscosity  (L2/f) 

Density  of  fluid  ( M/L 2) 

Shaft  angular  velocity  (1/t) 

Subscripts 

Reservoir  value,  radial  component 
Sump  value 
Tangential  component 


f5< 


Table  1.  Definition  of  symbols  used  in  figures. 


Supply  Pressure 

Rotor  Speeds 

Inlet  Circumferential  Velocities 

1  -  3.03  bar 

2  -  4.46  bar 

3  -  5.84  bar 

4  -  7.22  bar 

5  -  8.25  bar 

1  -  3,000  cpm 

2  -  6,000  cpm 

3  -  9,500  cpm 

4  -  13,000  cpm 

5  -  16,000  cpm 

1  -  Zero  circumferential  velocity 

2  -  Low  velocity  with  rotation 

3  -  High  velocity  with  rotation 

Honeycomb  Stator  Seals 


o'  o  o 
1  I  I 
-<f  m  <© 


Figure  3.  Seal  configurations. 
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Figure  9.  Cross-counted  atiffnesa  versus  inlet  circumferential  velocity  ratio  for  seals 
1-6  of  figure  3. 
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Figure  13.  A  comparison  of  experimental  and  theo-  Figure  14.  A  comparison  of  experimental  and  theo¬ 
retical  cross-coupled  stiffness  versus  ro-  retical  cross-coupled  stiffness  versus  in- 

tor  speed  for  seal  1  of  figure  3.  let  circumferential  ratio  for  seals  1-3  of 

figure  3. 


A  comparison  of  experimental  and  the-  Figure  16.  A  comparison  of  experimental  and  the¬ 
oretical  direct  damping  versus  rotor  oretical  direct  damping  versus  inlet  cir- 

speed  for  seal  1  of  figure  3.  cumferential  velocity  ratio  for  seals  1-3 

of  figure  3. 


THEORY  VERSUS  EXPERIMENT  POR  THE 
ROTORDYNAMIC  COEFFICIENTS  OF  LABYRINTH  GAS  8EALS: 
PART  I  -  A  TWO  CONTROL  VOLUME  MODEL1 


Joseph  K.  Schaerer 
Mechanical  Engineering  Department 
Texas  A&M  University 
College  Station,  Texas  77843 


SUMMARY 

The  b  asic  equations  are  derived  far  a  two-control- volume  model  for  compressible  flow  in  a  labyrinth  seal.  The 
recirculation  velocity  in  the  cavity  is  incorporated  into  the  model  for  the  first  time.  The  flow  is  assumed  to 
be  completely  turbulent  and  isoenergetic.  The  wall  friction  factors  are  determined  using  the  Blasius  formula. 
Jet  flow  theory  is  used  for  the  calculation  of  the  recirculation  velocity  in  the  cavity.  Linearised  seroth  and 
first-order  perturbation  equations  are  developed  for  small  motion  about  a  centered  position  by  an  expansion 
in  the  eccentricity  ratio.  The  seroth-order  pressure  distribution  is  found  by  satisfying  the  leakage  equation 
while  the  circumferential  velocity  distribution  is  determined  by  satisfying  the  momentum  equations.  The 
first-order  equations  are  solved  by  a  separation  of  variable  solution.  Integration  of  the  resultant  pressure 
distribution  along  and  around  the  seal  defines  the  reaction  force  developed  by  the  seal  and  the  corresponding 
dynamic  coefficients. 

INTRODUCTION 

The  problems  of  instability  and  synchronous  response  in  turbomachines  have  arisen  recently  because  of 
the  trends  in  design  toward  greater  efficiency  with  higher  performance.  To  achieve  these  design  goals,  the 
machine*  are  designed  for  higher  speeds,  larger  loadings,  and  tighter  clearances.  As  loadings  are  increased 
and  clearances  decreased,  fluid  forces  increase  and  can  lead  to  unstable  or  self-excited  vibrations.  One  of 
the  rotordynamk  force  mechanisms  which  has  been  shown  to  cause  self-excited  vibration  and  synchronous 
response  in  centrifugal  compressors  is  that  of  the  forces  developed  by  labyrinth  seals. 

The  flow  in  a  labyrinth  seal  has  been  shown  by  experiment  [l]  and  calculation  [2,3,4]  to  be  comprised  of  two 
flow  regimes:  a  jet  flow  region  in  the  leakage  path  and  a  recirculating  velocity  region  in  the  cavity  itself  (see 
figure  1).  The  first  attempts  at  analysis  of  this  system  neglected  the  axial  velocity  components  in  the  flow 
and  concentrated  on  the  circumferential  components.  This  was  the  single-control-volume  approach,  used  in 
refs  [5,6].  In  an  attempt  to  improve  upon  the  results  of  these  analyses,  the  two-control-volume  approach 
was  introduced.  The  most  notable  of  the  two-control-volume  analyses  is  that  of  Jenny  et  aL  [7].  Jenny  et 
aL  [7]  used  the  two-control-volume  approach  in  conjunction  with  a  2-D  CFD  solution  to  the  Navier-Stokes 
equations  to  account  for  the  free  shear  stress  between  the  jet  flow  and  the  cavity  flow.  However,  they 
neglected  the  recirculating  velocity  in  the  cavity  and  assumed  the  flow  to  be  incompressible.  Further,  the 
present  author  obtains  different  signs  in  the  expansion  of  the  continuity  equation  and  different  perturbation 
equations.  These  discrepancies  are  explained  in  detail  in  Appendix  D.  The  theory  of  Jenny  et  aL  [7]  showed 
substantial  improvement  in  the  prediction  of  stiffness  and  damping  coefficients,  but  in  the  end,  correction 
factors  had  to  be  incorporated  into  the  calculation  of  the  shear  stress  to  improve  the  correlation  with  test 
data. 

This  paper  introduces  the  calculation  of  the  recirculation  velocity  into  the  analysis.  The  model  for  the 
recirculation  velocity,  C/a,  used  here  is  illustrated  in  figure  2.  This  velocity  component  is  important  in  the 
calculation  of  the  cavity  shear  stresses.  The  focus  is  on  the  shear  stresses,  because  experimental  result  [8] 
have  shown  that  the  stiffness  and  damping  coefficients  are  very  sensitive  to  the  circumferential  velocity  in  the 
seaL  In  the  control  volume  analysis  to  be  presented,  the  solution  to  the  circumferential  momentum  equation 

1This  work  was  supported  in  part  by  NASA  Grant  NAS3-181  from  NASA  Lewis  Research  Center  (Technical 
Monitor,  Robert  Hendricks)  and  AFOSR  Contract  F49820-82-K-0033  (Technical  Monitor,  Tony  Amos) 


yields  the  circumferential  velocity  in  the  seal.  Ad  improvement  in  the  shear  stress  calculation  will  yield  an 
improvement  in  the  calculation  of  the  stiffness  and  damping  coefficients.  The  CFD  results  of  Rhode  |2,3}  are 
used  to  evaluate  the  models  for  shear  stress  and  recirculation  velocity  used  in  this  paper.  The  results  of  this 
analysis  are  compared  to  a  new  set  of  experimental  results  for  teeth-on-rotor  and  teeth-on-stator  labyrinth 
seals  in  Part  2  a ( this  paper. 

CONTROL  VOLUME  MODELLING 

Before  proceeding  with  the  solution  development,  the  approach  taken  in  modelling  the  flow  will  be  discussed. 
As  mentioned  previously,  the  flow  in  a  labyrinth  seal  is  known  to  have  two  distinct  regions:  a  jet  flow  region 
in  the  leakage  path  and  a  recirculating  flow  region  in  the  cavity  itself  (see  figure  l).  Therefore,  a  two-control- 
volume  model  seems  appropriate.  The  choice  of  control  volumes  is  between  the  ”box-in-a-boot*  model  (see 
figure  3)  of  Jenny  et  aL  (7)  or  a  more  conventional  model  with  a  control  volume  for  the  jet  flow  and  one 
for  the  recirculating  flow  in  the  cavity,  as  shown  in  figure  2.  The  two-separate-control- volume  model  was 
chosen  since  it  as  suggested  by  the  known  physics  of  the  flow.  The  flow  enters  the  seal  and  separates  into 
two  distinct  flow  regions  which  are  separated  by  the  dividing  streamline. 

The  final  question  is  whether  the  control  volumes  should  be  defined  using  a  geometric  boundary  or  using 
the  dividing  streamline  as  the  boundary.  The  dividing  streamline  approach  seems,  at  first,  to  be  the  obvious 
choice.  The  governing  equations  would  be  simplified  by  the  restriction  of  no  flow  across  a  streamline,  the 
free  shear  stress  relations  are  derived  for  flow  along  the  dividing  streamline,  and  the  solution  for  the  velocity 
of  the  recirculating  flow  may  be  derived  for  flow  along  the  dividing  streamline.  Despite  these  advantages,  the 
dividing  streamline  approach  was  not  used,  due  to  mathematical  constraints.  The  mathematical  limitations 
of  the  dividing  streamline  approach  are  dealt  with  in  detail  by  Scharrer  (0].  The  geometric  boundary  approach 
relies  on  the  assumption  that  the  dividing  streamline  makes  a  small  angle  with  the  horisontaL  As  will  be 
shown,  this  is  a  good  assumption  which  has  been  verified  experimentally.  The  geometric  boundary  approach 
and  solution  is  provided  in  the  following  section. 

GEOMETRIC  BOUNDARY  APPROACH 
Aaaumptions 

(1)  The  fluid  is  considered  to  be  an  ideal  gas. 

(2)  Pressure  variations  within  a  chamber  are  small  compared  to  the  pressure  differences  across  a  seal 
strip. 

(3)  The  lowest  frequency  of  acoustic  resonance  in  the  cavity  is  much  higher  that  that  of  the  rotor  speed. 

(4)  The  eccentricity  of  the  rotor  is  small  compared  to  the  radial  seal  clearance. 

(5)  Although  the  shear  stress  is  significant  in  the  determination  of  the  flow  parameters  (velocity  etc.  ), 
the  contribution  of  the  shear  stress  to  the  forces  on  the  rotor  are  negligible  when  compared  to  the  pressure 
farces. 

(6)  The  cavity  flow  is  turbulent  and  isoenergetic. 

(7)  The  recirculation  velocity,  (7a,  is  unchanged  as  it  swirls  within  a  cavity. 

Procedure 

The  analysis  presented  here  is  developed  for  the  teeth-on-rotor  *  see-through”  labyrinth  seal  shown  in  figure  5. 
The  equivalent  equations  for  the  teeth-on-stator  labyrinth  seal  are  given  in  Appendix  A.  The  continuity  and 
circumferential  momentum  equations  will  be  derived  for  the  two-control-volume  model  shown  in  figures  2,6,7 
and  8.  A  leakage  model  will  be  employed  to  account  for  the  axial  flow.  The  governing  equations  are  linearised 
using  a  perturbation  analysis  for  small  motion  about  a  centered  position.  The  seroth-order  continuity  and 
momentum  equations  will  be  solved  to  determine  the  steady  state  pressure,  axial  and  circumferential  velocity 
for  each  cavity.  The  first-order  continuity  and  momentum  equations  will  be  reduced  to  linearly  independent, 
algebraic  equations  by  assuming  an  elliptical  orbit  for  the  shaft  and  a  corresponding  harmonic  response  for 
the  pressure  and  velocity  perturbations.  The  force  coefficients  for  the  seal  are  found  by  integration  of  the 
first-order  pressure  perturbation  along  and  around  the  shaft. 

GOVERNING  EQUATIONS 
Continuity  Equation 

The  control  volumes  of  figures  2  and  6  have  a  unity  circumferential  width.  Their  continuity  equations  are: 


dpA\  dpWxAi 

"ir+  Rt~dT+nK+t~mi+mr=0  w 

BpA2  8pW2A2  . 

dt  +  R»2ae  "v  0  (2> 

For  the  teeth-on-rotor  cue,  A\  =LCr,  A2  =  LB,  /?«i  —  Rs,  and  R»2  -  Re  +  B. 

Momentum  Equation* 

The  following  momentum  equation*  for  control  volume*  I  and  II  are  derived  naing  figure*  7  and  8  which 
•how  the  pressure  force*  and  ahear  atreaae*  acting  on  the  control  volume*. 

9pWlAl  2pWiA\  8WX  pWx  8Ax  ,  WXAX  dp  .  . 
at  +  Rtl  86  Rex  86  +  Rex  86  + 

A  dP' 

+  m+iWxi  -  fh.Wi-1  =  (3) 


8PW2A2  lpW2A28W2  pW2  8A2  W2A28p 
dt  +  Re 3  86  +  Re2  86  +  Re2  86 

+  rr^Wd.  -  ~  rj*Li  +  V» ari L>i  (4) 

where  ar  and  a*  are  the  dimenaionleaa  length  upon  which  the  shear  atreaae*  act  and  are  defined  for  the 
teeth-on-rotor  labyrinth  by 

«*,  =1  ar{  -  (2 Bi  +  Li)/Li  (5) 

We ,•  la  the  circumferential  velocity  between  the  control  volumes. 

Blaaius  [10]  determined  that  the  shear  stresses  for  turbulent  flow  in  a  smooth  pipe  could  be  written  as 


where  is  the  mean  flow  velocity  relative  to  the  surface  upon  which  the  shear  stress  is  acting.  The  constants 
mo  and  no  can  be  empirically  determined  for  a  given  surface  from  pressure  flow  experiments.  However,  for 
smooth  surfaces  the  coefficient*  given  by  Yam  ad*  [ll]  for  turbulent  flow  between  mimbr  surfaces  are: 

mo  *=  —0.25  no  =  0.079 

Applying  Blaaius’  equation  to  the  labyrinth  rotor  surfaces  yields  the  following  definitions  for  the  rotor  shear 
stress  in  the  circumferential  direction.  Note  that  the  recirculation  velocity,  U2,  is  included  in  the  definition 
of  the  total  velocity  acting  on  the  rotor. 


rr  =  2  ~  W*)7  +  ~  W») nr  I 


(Rejv  -  W3 


where  D^2i  is  the  hydraulic  diamter  of  control  volume  II,  defined  by 

Dhv  —  2  BiLi/{Bi  +  Li) 

Similarly,  the  stator  shear  stress  in  the  circumferential  direction  is: 


y/wl  +  U?Dh 


where  D^u  is  the  hydraulic  diameter  of  control  volume  I,  defined  by 

Dku  =  2CV,L,/(Cr.  +  A)  (9) 

and  the  axial  velocity,  Ult  is 

U\  =  m  j  pC  r,  (10) 

Figure  9  shows  a  comparison  of  the  predictions  from  equation  (8)  and  CFD  results  for  stator  wall  shear  stress 
for  seal  A  of  table  1.  The  figure  shows  that  the  comparison  is  very  good.  Similar  results  are  obtained  for  the 
other  seals  of  table  1.  Figure  10  shows  a  comparison  of  rotor  wall-shear-stress  predictions  from  equation  (8), 
CFD  and  averaged  CFD  for  seal  A  of  table  1.  The  averaged  CFD  results  is  used  here  for  comparison  since 
the  bulk  flow  model  yields  a  single  averaged  result  for  cavity  shear  stress  and  is  not  capable  of  modelling 
the  complex  flowfield.The  figure  shows  that  the  prediction  of  equation  (6)  is  close  to  the  CFD  results.The 
dips  in  the  CFD  results  are  the  lower  corners  of  the  cavity.  Similar  results  are  obtained  for  the  other  seals 
of  table  1. 

Table  1.  Seal  geometries  calculated  by  Rhode. 

Seal 


A 

B 

C 

D 

Rs 

72.5043mm 

72.5043mm 

72.5043mm 

41.780mm 

B 

3.175mm 

S.175mm 

S.175mm 

0.889mm 

L 

3.175mm 

S.175mm 

3.969mm 

0.8585mm 

TP 

0.35mm 

0.35mm 

0.35mm 

0.15mm 

Cr 

0.4064mm 

0.508mm 

0.508mm 

0.2159mm 

The  flow  across  a  labyrinth  tooth  is  very  similar  to  the  flow  of  a  turbulent  jet  issuing  from  a  slot.  The  problem 
with  using  jet-flow  results  for  labyrinth  seals  is  that  current  jet-flow  theory  only  considers  the  flow  of  a  jet  with 
a  coflowing  stream  or  a  crossflowing  stream,  not  both.  In  the  following  derivation,  the  relationships  given  by 
Abramovich  {12]  for  the  velocity  profile  of  a  semi-contained,  one-dimensional,  turbulent  jet  with  a  coflowing 
stream  are  assumed  to  apply  for  the  two-dimensional  labyrinth  seal  flow.  According  to  Abramovich  (12],  the 
velocity  profile  for  such  a  flow  can  be  shown  to  fit  the  following  function  when  compared  to  experimental 
results: 

v  =  +  (vj  -  v,)  |i  -  |  (11) 

where  the  coordinate  y,  the  mixing  thickness  b,  and  the  boundary  layer  thickness  are  defined  in  figure 
11.  The  relationship  between  the  boundary  layer  thickness  and  the  mixing  thickness  was  found  (12]  by 
comparison  to  experiment  to  be: 

ya/fc  =  0.584  -  0.134(wa/wi)  (12) 

Once  the  velocity  ratio  across  the  dividing  streamline,  vj/vi,  is  found,  equation  (12)  reduces  to  a  constant. 
The  total  free  shear  stress  is  found  using  Prandtl’s  mixing  length  hypothesis  (13]: 

'/«  =  P? 

where  the  mixing  length,  t,  for  a  labyrinth  seal,  has  been  determined  from  the  calculations  of  Rhode  (2,3]  to 
be: 

t  =  0.275 b  (14) 

Table  1  shows  the  seal  geometries  calculated  by  Rhode  (2,3].  The  mixing  length,  l,  given  in  equation  (14) 
is  the  most  sensitive  factor  in  this  solution.  The  large  magnitude  of  the  mixing  length  shows  the  high 


turbulence  level  of  the  labyrinth  flow  ae  compared  to  aimilar  flow*.  The  typical  values  given  for  the  mixing 
lengths  of  rectangular  and  round  jet  flows,  in  one  dimension,  are  in  the  range  of  0.07  to  0.09.  Without  the 
CFD  results,  one  of  these  values  would  have  to  be  used  and  the  results  of  using  t/b  in  the  range  [0.07,0.09j 
would  have  been  disappointing. 

Jenny  et  al.  (7]  used  a  2-D  CFD  code  to  obtain  a  correlation  for  t/b  as  a  function  of  clearance  and  tooth 
geometry.  Their  relation  is  shown  below  for  the  teeth-on-rotor  case: 

t/b  =  0.055(1  ■+  1.030/1  +  0.08  \/Rs/L)  (15) 

However,  their  shear  stress  relation  neglected  the  recirculating  velocity  component,  U2.  Upon  comparison 
with  the  data  of  Rhode  (2,3],  the  mixing  length  ratio,  t/b,  was  found  to  be  relatively  constant  when  the 
shear  stress  is  calculated  using  aD  velocity  components. 

Substituting  the  differentiated  version  of  equation  (11)  and  equation  (14)  into  equation  (13)  yields  an  ex¬ 
pression  for  the  total  free  shear  stress.  At  the  interface  of  the  two  control  volumes  (y=0),  the  total  free  shear 
stress  is: 

rit  *  0.68,1(03  -  Vl)l  [l  -  (^) l'6]  3  (^)  (16) 

The  circumferential  component  of  the  free  shear  stress  is: 

iy  =  0.68,v/(W3-Wa)>  +  (I/3-^)a  [W2  -  W,)  [l  -  (^)  *  (^)  (17) 

The  circumferential  component  at  velocity  at  the  interface,  W^,  is  obtained  from  equation  (11). 


v.  +  <*,-v,)[i -(£)*•*]’ (E) 


Equations  (16,17,18)  are  all  valid  along  the  dividing  streamline.  Since  the  control  volumes  are  defined 
geometrically  and  not  by  the  dividing  streamline,  the  shear  stress  calculated  using  the  above  equations  is 
assumed  to  be  close  to  that  existing  along  the  geometric  boundary  line.  This  is  a  good  assumption  considering 
that  the  angle  of  the  dividing  streamline  from  the  horisontal  has  been  found  experimentally  to  be  on  the 
order  of  6  degrees  by  several  investigators  (14,15], 

The  analysis  to  this  point  is  incomplete  in  that  the  recirculation  velocity  ,  U2,  and  the  relationship  between 
the  mixing  thirhness  and  the  boundary  layer  thickness,  Va/b,  are  undefined.  The  following  section  presents 
the  jet  flow  theory  used  to  determine  the  recirculation  velocity,  U2,  and  subsequently  yta/b.  The  discussion 
is  rather  lengthy,  but  the  final  result  is  relatively  simple. 

Determination  of  tke  rtcirc ulatio a  velocity 

As  mentioned  previously,  there  is  no  closed-form  solution  for  the  flowfield  present  in  a  labyrinth  seal  cavity. 
The  flow  is  highly  three  dimensional  and  completely  turbulent.  An  approximation  for  the  velocity  profile 
can  be  obtained  using  the  theory  for  the  flow  of  a  two-dimensional,  turbulent,  isoenergetic,  half-infinite  jet. 
Figure  12  shows  the  model  for  this  theory.  The  flow  is  assumed  to  enter  with  one  velocity  component,  in  the 
x-direction,  and  spread  into  the  cavity  developing  a  y-component  of  velocity.  The  model  doe#  not  account 
for  the  circumferential  velocity  component  which,  is  the  same  order  of  magnitude  as  the  »*ial  velocity,  in  a 
labyrinth  seal  flowfield.  The  solution  procedure  involves  solving  the  infioHeasimal  form  of  the  x-momentum 
equation  for  the  dimensionless  velocity  profile  and  then  solving  the  integral  form  of  the  continuity  and 
momentum  equations  in  order  to  determine  the  dimensionless  velocity  along  the  dividing  streamline. 

The  following  is  a  summary  of  the  derivation  of  the  equations  necessary  to  determine  the  dimensionless 
velocity  along  the  dividing  streamline.  A  complete  discussion  of  this  theory  can  be  found  in  Korst  et  al.  (I6j 
and  Scharrer  (9].  The  following  derivation  uses  the  assumption  that  the  curvature  in  the  dividing  streamline 
is  small  The  infinitesaimal  form  of  the  x-momentum  equation  which  has  been  reduced  using  the  continuity 
equation  is: 


(19) 


du  da  d  f  du\ 

fur,+f’s;‘a;{‘l,ai) 


where  c  ia  the  apparent  (turbulent)  kinematic  viscosity  and  the  x  and  y  velocity  components,  u  and  v, 
respectively,  are  time  averaged.  Since  the  flow  illustrated  in  figure  12  is  a  quasi-one-dimensional  jet  flow 
where  there  is  little  or  no  initial  vertical  velocity  component,  equation  (19)  can  be  linearised  using  a  small 
perturbation  method.  The  following  simplified  equation  of  motion  is  obtained: 

„  da"  a3u" 

Ul  ds  ~  *  dy1  ^ 

Here  the  turbulent  viscosity  is  expressed  in  a  modified  form  of  Prandtl’s  exchange  coefficient,  e,  so  that  after 
introducing  the  dimensionless  variables 

9  V .  Um 

f-f 

t  (21) 

*  =  7 


+6Uifl+) 
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where 


1.0  as 


and  by  the  transformation 


e_  f*  +/{*)*+ 
*  Jo  2o* 

the  following  formulation  of  the  equation  of  motion  is  obtained 


with  the  initial  conditions 


a 4>  d2* 
at  ~  a fa 


4  =  ^(Q|  f)  =  o 
4  =  4(o,<)  =  44<) 
4  =  ^(o,  f)  =  l.o 


for  —  oo  <  f  <  0 
for  0  <  f  <  1.0 
for  1.0  <  f  <  oo 


and  boundary  conditions 


4  =  4{Z,~°°)  ~  0  for  (  <0 
4  =  4(t,  oo)=1.0  for  £  >  0 

Here  6  is  the  initial  boundary  layer  thickness  shown  in  figure  12.  o  is  the  jet  spreading  parameter  which  has 
been  found  experimentally  by  Korst  and  IHpp  (17)  to  fit  the  following  equation: 

o  =  12.0  -f  2.758Af#  for  air  (23) 

The  solution  to  equation  (22)  for  the  above  initial  and  boundary  conditions  is: 


4  =  0.5(1  -er/(*-,)J+4-r 

V*  Jn-fi, 


where  rjp,  the  position  parameter  ia  given  by 


Vp  =  ^  and  1  =  (V  p 
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The  error  fraction,  crf(x),  i>  defined  by 


where  erf(-x)  =  —erf(x) 


Looking  at  *  limiting  cue  of  equation  (24): 


z 

at  —  — *  oo 
0 


then  ip  — — *  oo 
{  — *  oo 


Thin  limiting  cate  it  for  either  no  initial  boundary  layer,  which  it  a  good  aaeumption  for  labyrinth  teals,  or 
folly  developed  velocity  profile*.  Since  qp  — *  0,  the  variable  fj  it  now  ondefined.  Liepman  and  Laofer  (18] 
have  shown  that  rj  for  thit  limiting  condition  it: 


(25) 


Goertler  (19]  hat  eh  own  that  the  dimenaionlett  velocity  ,4>,  follows  directly  from  equation  (24)  when  r)p  — »  0: 


*  =  0.5(1  +  erf{V))  (26) 

Equation  (26)  it  a  solution  for  the  dimensionless  velocity  profile,  p,  at  any  dimensionless  position,!).  The 
goal  of  thi*  development  it  to  determine  the  dimensionless  velocity,  4>i,  along  the  dividing  streamline  whose 
dimensionless  position  is  fy.  The  dimensionless  dividing  streamline  position,  tj}  ,  can  be  obtained  by  solving 
the  integral  form  of  the  continuity  and  x-momentum  equations  for  the  system  shown  in  figure  12.  Equations 
(23,25  and  26)  are  used  to  obtain  the  solution  to  the  integral  equations  which  are  to  be  derived  next. 

Control  Volume  Analytu 

The  coordinate  systems  and  definition  of  the  control  surface  are  shown  in  figure  12.  The  (x,y)  coordinate 
system  is  the  intrinsic  coordinate  system  while  the  (X,Y)  coordinate  system  is  the  reference  system. Equation* 
(25)  and  (26)  are  approximate  relations;  exact  relationships,  if  known,  would  provide  for  conservation  of 
momentum  for  the  constant  pressure  mixing  region.  The  reference  coordinate  system  is  the  coordinate 
system  in  which  momentum  is  conserved.  The  intrinsic  coordinate  system  is  located  with  respect  to  the 
reference  coordinate  system  by  a  control  volume  analysis  utilising  the  conservation  of  momentum  principle 
for  this  constant  pressure  mixing  region.  The  relationship  between  the  coordinate  systems  normal  to  the  jet 


¥«(*)  =  V~Y  vith  Vm(0)  =  0. 

X -Momentum  Equation 

The  steady  flow  x-momentum  equation  for  the  jet  flow  shown  in  figure  12,  written  for  the  intrinsic  coordinate 
system  and  expressed  in  the  previously  defined  dimensioless  variables  is: 

V*  f  +  f/R  “  Vp  =  f  —i7drf  +  i)m  (27) 

JO  P l  J- oo  Pi 

where  distance  R  is  chosen  such  that: 

l  -4>{vr)  «<  i.o 

Far  the  momentum  equation,  the  lower  control  surface  is  located  at  —  oo.  This  equation  contains  no  surface 
forces.  This  is  realistic  for  a  labyrinth  seal  if  location  R  is  chosen  far  from  the  stator  walL  Applying  the 
condition  of  no  initial  boundary  condition  (t)p  — *  0)  equation  (27)  is: 


[I 

[4 


Continuity  Equation 

The  steady  flow  continuity  equation,  written  for  the  intrinsic  coordinate  system  and  expressed  in  the  previ¬ 
ously  defined  dimensionless  variables  is: 


f10  p  p 

/  —4>o<k  +  fJR-flp-  /  —4>dr\  +  r?m 

J  0  Pi  Jo,  Pi 


For  the  continuity  equation, the  lower  control  surface  is  coincident  with  the  jet  dividing  streamline.  Substi¬ 
tuting  the  results  of  the  momentum  equation,  equation  (27),  into  equation  (28)  yields: 


r  -*dr, = [l0  £(1  -  *0  )ds  +  r  U'i* 

Jn,  Pl  Jo  Pi  J-oa  Pi 


Making  the  assumption  of  no  initial  boundary  layer  (qp  — »  0),  equation  (SO)  becomes: 


r  =  r 

Jm,  Pl  J- OO  Pl 


The  density  ratio,  (p/ pi),  for  isoenergetic  flow  (constant  temperature)  is  given 

P  (1  ~  Ca7) 
pi  (1  -  Ca7*7) 

The  final  form  of  the  continuity  equation  becomes: 


/"•*  *  t  f**  *7 

Jn,  (l-CaV*)*'";—  (1-CoV3)^ 


where  Ca  is  the  Crocco  number.  liquation  (S3)  can  be  numerically  integrated  for  a  given  Croce o  number 
using  the  definition  for  *  given  in  equation  (26).  The  Crocco  number  is  defined  as: 


Ca?  = 


(7  —  i)M* 

(2  +  (7-l)AP) 


The  Crocco  number  is  a  dimensionless  velocity  similar  to  the  Mach  number.  The  Crocco  number  uses  the 
maximum  isen tropic  speed  of  a  gas  while  the  Mach  number  uses  the  local  speed  of  sound.  The  Mach  number 
varies  between  0  andoo  while  the  Crocco  number  has  a  range  of  0  to  1. 

The  solution  to  equation  (S3), the  location  of  the  dividing  streamline,  can  be  obtained  by  the  foDowing  steps: 

0)  Calculate  the  Mach  number  using  the  seroth-order  leakage  value. The  seroth-order  leakage  is  discussed 
in  the  next  section. 

1)  Calculate  the  Crocco  number  using  equation  (34). 

2)  Substitute  equation  (26)  into  equation  (33)  and  integrate  the  error  function.  The  value  of  the  error 
function  at  the  limits  R  and  —00  is  1.0,  leaving  an  equation  in  9/  only.  This  is  solved  for  9/  which  is  the 
dimensionless  location  of  the  dividing  streamline. 

3)  Insert  rj}  into  equation  (26)  to  obtain  the  dimensionless  velocity  along  the  dividing  streamline, 

The  results  of  this  solution  procedure  are  tabulated  b  table  2,  for  air.  For  air  (7  •=  L4)  flowing  in  a  labyrinth 
seal,  the  maximum  possible  Mach  number  is  1.0.  Therefore,  the  maximum  possible  Crocco  number  is  0.408 
or  Ca 3  =  0.167.  The  range  of  solutions  is: 

0.61632  <  dy  <  0.6263 

Using  an  average  solution  of  =  0.62  gives  a  maximum  error  of  less  th  n  ±1  percent.  The  recirculation 
velocity  at  the  interface  is: 


i 


t 


a 

A 


'A 


U2j  =  0.«2t/i  (35) 

The  only  remaining  problem  is  the  numerical  definition  of  yj/6.  Looking  back,  equation  (11)  and  equation 
(26)  both  describe  the  axial  velocity  profile  in  the  jet  flow  field.  If  the  following  observation  is  made 


then  eqnation  (35)  can  be  substituted  back  into  equation  (12)  yielding  the  following  numerical  definition  for 
Vi/b: 


=  0.564  -  0.134*,  «=  0.5 

V 

It  is  interesting  to  note  that  Jenny  et  aL  (7j  assumed  that  y^/b  —  0.5. 

Figure  13  shows  a  plot  of  the  dimensionless  axial  velocity  profile  in  the  recirculation  region  for  seal  A  of 
table  1  as  calculated  by  Rhode  (2,3].  This  profile  is  for  the  center  of  the  recirculation  region  to  the  top  of 
the  labyrinth  tooth.  The  intersection  of  the  two  dashed  lines  is  the  location  and  value  of  the  theoretical 
recirculation  velocity  as  calculated  using  equation  (35)  and  the  assumption  that  the  dividing  streamline 
makes  an  angle  of  <1°  with  the  horisontaL  The  agreement  is  excellent.  Again,  equation  (35)  was  derived 
using  a  2-D  theory  which  neglects  the  circumferential  velocity  component.  Equation  (35)  is  actually  the 
velocity  at  the  interface  of  the  two  control  volumes.  The  velocity  components  used  ia  the  shear  stress 
equations  are  all  average  velocity  components.  Tb  be  consistent,  the  average  recirculation  velocity  must  be 
used.  The  CFD  results  show  that  the  velocity  distribution  is  parabolic  in  nature.  Integrating  this  yields: 


U2  =  0.206 J/j  (36) 


Table  2.  Tabulated  solution  to 

equation  (33). 

Ca 7 

Co 7 

0.00000 

0.61632 

0.68000 

0.67553 

0.05000 

0.61815 

0.72000 

0.68188 

0.10000 

0.62211 

0.76000 

0.68803 

0.15000 

0.62523 

0.80000 

0.68724 

0.20000 

0.62848 

0.84000 

0.70688 

0.24000 

0.63128 

0.86480 

0.713844 

0.28000 

0.63405 

0.88360 

0.718844 

0.32000 

0.63725 

0.80250 

0.726834 

0.36000 

0.64047 

0.82160 

0.734848 

0.40000 

0.64387 

0.84080 

0.744883 

0.44000 

0.64748 

0.86040 

0.757869 

0.48000 

0.65132 

0.88010 

0.777432 

0.52000 

0.65543 

0.882016 

0.788766 

0.56000 

0.65878 

0.888001 

0.823427 

0.60000 

0.66462 

1.000000 

1.000000 

0.64000 

0.66882 

Reduced  Equation* 

The  solution  of  the  governing  equations  can  be  simplified  by  reducing  the  number  of  equations  by  one.  This 
reduction  is  accomplished  by  using  equation  (2)  to  eliminate  mr  from  the  other  equations.  The  continuity 
equation  for  control  volume  I  becomes: 


BpA\  t  dpWiAi  _  ,  dpA2  dW7A2  n 

+  +  +: 1 "  mi +  ~dT  +  rm7  =  0 


(37) 


If  equation  (37)  times  the  circumferential  velocity,  Wx  is  now  subtracted  from  equation  (3),  the  following 
reduced  form  of  the  momentum  equation  for  control  volume  I  is  obtained: 


y.V. 


pWxAidWi  dpAi  .  8pW7A7 


+  -  Wi.-i)  =  +  *j.£.  - 


Similarly,  if  equation  (2)  time*  the  circumferential  velocity,  W7,  ie  eubtracted  from  equation  (4),  the  reduced 
momentum  equation  for  control  volume  II  ii  obtained. 

aw7  pW7A7aw7  (dpA7  dpW7A7 


=  ~  Rm7  88  +T,<Li  TriaTiLi 

The  number  of  variable  ia  reduced  by  using  the  ideal  gas  law  to  eliminate  the  density  terms. 

Pi  =  PiRT  (40) 

Leakage  Equation 

To  account  for  the  leakage  mass  flow  rate  in  the  continuity  and  momentum  equations,  the  following  model 
was  chosen. 


m,  =  puPvHi\ 


lP?-r~P? 


where  the  kinetic  energy  carryover  coefficient,  pj  is  defined  by  Vermes  [20]  for  straight  through  scab  as: 

<■ 


where  a  — 


(Ihged  +  w) 


and  is  unity,  by  definition,  for  the  first  tooth  of  any  seal  and  all  the  teeth  in  interlocking  and  combination 
groove  seals.  This  definition  of  the  carryover  coefficient  is  a  local  coefficient  which  can  be  perturbed  in  the 
clearance.  The  previous  analyses  by  Childs  and  Schaner  [6]  and  Jenny  et  al  [7]  used  a  global  definition 
which  could  not  be  perturbed. 

The  flow  coefficient  is  defined  by  Chaplygin  [21]  as: 


A 

A 


i 


v. 


a. 
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where  «  =  «„/Cr,  is  the  eccentricity  ratio.  The  seroth-order  equation*  define  the  leakage  man  flow  rate  and 
the  circumferential  velocity  distribution  for  a  centered  position.  The  first-order  equations  define  the  per¬ 
turbations  in  pressor*  and  circumferential  velocity  due  to  radial  position  perturbation  of  the  rotor.  Strictly 
speaking,  results  at  a  first  order  analysis  are  only  valid  for  small  motion  about  a  ceutered  position. 

Zcroth-Order  Solution 

The  seroth-order  leakage  equation  is 


rrv,+1  =  m,  +  m0  (45) 

and  is  used  to  determine  both  the  leakage-rate,  nv,  and  pressure  distribution  for  a  centered  position.  The 
leakage-rate  and  cavity  pressures  are  determined  iteratively,  in  the  following  manner.  First,  determine 
whether  the  flow  is  choked  or  not  by  assuming  that  the  Mach  number  at  the  last  tooth  is  one.  Then,  knowing 
the  pressure  ratio  for  flow  at  sonic  conditions,  the  pressure  in  the  last  cavity  is  found.  The  mass  flow  can 
be  calculated  using  equation  (44).  Working  backwards  towards  the  first  tooth,  the  rest  of  the  pressures  can 
be  found  using  equation  (41).  The  final  pressure  calculation  will  result  in  the  reservoir  pressure  necessary 
to  produce  the  sonic  condition  at  the  last  tooth.  If  the  actual  reservoir  pressure  is  less  than  this  value,  then 
the  flow  is  unchoked.  Otherwise,  it  ia  choked.  If  the  flow  is  choked,  a  similar  procedure  is  followed,  but  now 
the  pressure  in  the  last  cavity  is  guessed  and  a  mass  flow  rate  calculated  using  equation  (44).  The  remaining 
pressures  are  calculated  using  equation  (41).  This  is  repeated  until  the  calculated  reservoir  pressure  equals 
the  actual  reservoir  pressure.  If  the  flow  is  unchoked,  the  pressure  in  the  first  cavity  is  guessed  and  a  mass 
flow  rate  calculated  using  equation  (41).  The  remaining  pressures  are  calculated  with  the  same  equation. 
This  procedure  is  repeated  until  the  calculated  sump  pressure  equals  the  actual  sump  pressure. 

The  seroth-order  circumferential-momentum  equations  are 


m«,(Wloi  -  Wl0i.l)  =  (r/oj  -  r.o.a*,)L, 


(46) 

rjoili  -  rroiariLi  (47) 

From  calculated  pressures,  the  densities  can  be  calculated  at  each  cavity  from  equation  (46),  and  the  only 
unknowns  remaining  in  equations  (46)  and  (47)  are  the  circumferential  velocities  Wloi  and  Given  an 
inlet  tangential  velocity,  a  Newton-root-finding  approach  can  be  used  to  solve  equations  (46)  and  (47)  for 
the  i-th  velocities,  one  cavity  at  a  time.  This  is  done  starting  at  the  first  cavity  and  working  downstream. 

Pint-Order  Solution 

The  governing  first-order  equations  (48,49,50),  define  the  pressure  and  velocity  fluctuations  resulting  from 
the  seal  clearance  function.  The  continuity  and  momentum  equations  follow  in  order: 


y 

Gu  dt  +<7*  89  +G*  86  +G*  86  +  GbiPli 

+  GniPn-i  +  GjiPu+i  =  —GtiHu  —  G&  —  Gux—gjf- 

(48) 

V, 

v  8Wiu  ,  XuWloi  8Wui  ,  [v  ,  X*W3o,‘|  8Pu  ,  v  8Pu 

Xu  dt  1  Rs,  86  '  l**  '  R.3  J  86  +X *  dt 

V, 

+  9-~L  +  *4.  A,  +  XuPu-x  +  X*Wui  +  Xr<W3U 

V 

—  rhoiWjK.x  «=  XtiHu 

(49) 

V 

„  dW3U  .  \Y*P*  .  YuW^]  dW3U  .  L  .  IfcWwl  8PU 

Yu  dt  f  [  JU3  '  Rs3  J  86  '  [r*  1  R,3  j  86 

>* 

+  n,Pi.  +  YuW3U  +  YmPu- 1  +  YnWlu  =  YtiHu 

ot 

(50) 

where  the  X’s,  Y’s  and  G’s  are  defined  in  Appendix  B.  These  perturbation  equations  are  very  different  from 
those  of  Jenny  et  al.  (7|  because  their  analysis  neglects  pressure  perturbations  in  the  leakage  and  shear  stress 
equations,  and  assumes  that  the  density  is  constant. 


►WWTWI  Twni 


If  the  shift  center  move*  in  in  elliptical  orbit,  then  the  Mil  clearance  function  cm  be  defined  aa: 


tH j  =  —a  cotwt  cotff  —  btirxwt  rinff 

Q  fa 

=  — — |co«{tf  —  wt)  +  eos(0  +  wt)]  —  —\eo$(8  —  wt)  -  co*(8  ■+■  wt)] 

*  * 

The  pressure  and  velocity  fluctuations  can  now  be  stated  in  the  associated  solution  format: 


Pu  —  P£  eot(6  +  wt)  +  P%  »\n(8  +  wt)  +  P*  «os(tf  —  wt)  •+■  P~  rin(8  —  wt)  (52) 

1VU,  =  W^eico»{8  +  wt)  +  W^st'f»(tf  +  wt)  +  W^eieo$(8  —  wt)  +  W^stnftf  -  wt)  (53) 

Win  =  W£jCot(8  +  wt)  -f  W^,r»n(0  +  wt)  +  Wj ^eos(0  —  wt)  +  W^»in(8  —  wt)  (54) 

Substituting  equations  (51),  (52),  (53)  and  (54)  into  equations  (43),  (49)  and  (60)  and  grouping  like  terms 
of  sines  and  cosines  (as  shown  in  Appendix  C)  eliminates  the  time  and  theta  dependency  and  yields  twelve 
linear  algebraic  equations  per  cavity.  The  resulting  system  of  equations  for  the  i-th  cavity  can  be  stated: 

+  [*](*.)  +  M*m1(X*i)  -  f  (ifc)  +  ;(Ci)  (55) 

where 

(x<)  «  {p+,pZ,K,p~,  w&t  wfc.  w~{,  w^,  w&,  w&,  w^, 

(x.+1)  =  {p&xtrt+»p«+fKi+x&i«+»wiu+x  *  ^ioi+n  ^3m'+1i  Wm+1«  **2«'+1i  ^2e.+l)T 

The  A  matrices  and  column  vectors  B  and  C  are  given  in  Appendix  C.  To  use  equation  (55)  for  the  entire 
solution,  a  system  matrix  can  be  formed  which  is  block  tridiagonal  in  the  A  matrices.  The  sise  of  this 
resultant  matrix  is  (12NC  X  12NC)  since  pressure  and  velocity  perturbations  at  the  inlet  and  the  exit  are 
assumed  to  be  sero.  This  system  is  easily  solved  by  various  linear  equation  algorithms,  and  yields  a  solution 
of  the  form¬ 
's-  + 

P*  -  ~FZi  +  j*** 

a  b  W 

P£- + 

K  «  f  +  7^ 

DETERMINATION  OF  DYNAMIC  COEFFICIENT 

The  force-motion  equations  for  a  labyrinth  seal  are  assumed  to  be  of  the  form: 


-Ccl-I-*  i](f)*[S  :]{?) 


The  solution  of  equation  (57)  for  the  stiffness  and  damping  coefficients  is  the  objective  of  the  current  analysis. 
The  solution  procedure  used  for  this  analysis  is  the  same  one  used  by  ChOds  and  Scharrer  [6j.  The  desired 
solution  for  the  stiffness  and  damping  coefficients  is: 


;%c  (s«) 

»  nc 

SOLUTION  PROCEDURE  SUMMARY 

In  review,  the  eolation  procedure  aeee  the  following  sequential  steps: 

a)  Determination  of  whether  the  flow  is  choked  or  not  using  equations  (41)  and  (44). 

b)  The  steady-state  pressure  distribution  and  leakage  are  found  using  equation  (41)  and/or  (44). 

c)  The  steady-state  circumferential  velocity  distribution  is  determined  using  equations  (40)  and  (47). 

d)  A  system  equation  is  formed  for  the  first-order  perturbation  variables  and  solved  using  the  cavity 
equation  (&5). 

e)  Results  of  this  first-order  perturbation  solution,  as  defined  in  equations  (56),  are  inserted  into  equation 
(58)  to  defined  the  rotordynamic  coefficients. 

CONCLUSIONS 

This  paper  has  presented  a  new  two-control-volume  analysis  for  the  rotordynamic  coefficients  of  labyrinth 
gas  seals  which,  for  the  first  time,  accounts  far  the  recirculation  velocity  in  the  seal  cavity.  The  analysis 
was  developed  in  conjunction  with  a  2-D  CFD  model  which  was  used  to  verify  the  shear  stress  and  jet  flow 
models  used.  A  comparison  between  the  CFD  results  and  the  results  from  the  “bulk  flow”  model  of  this 
paper  showed  the  following: 

1)  The  new  two-control-volume  model  accurately  predicts  the  stator  wall  shear  stress  for  a  teeth-on-rotor 
labyrinth  seal  cavity. 

2)  The  new  model  predicts  the  cavity  wall  shear  stress  within  25  percent  of  the  CFD  results  for  a 
teeth-on-rotor  labyrinth  seal. 

S)  The  2-D  jet  flow  theory  used  in  the  new  model  accurately  predicts  the  magnitude  of  the  recirculation 
velocity  along  the  dividing  streamline. 

4)  The  CFD  results  show  that  the  mixing  length  parameter,  L,  used  in  the  free  shear  stress  equation  is 
relatively  constant  and  need  not  be  considered  a  function  of  cavity  geometry  as  was  assumed  by  Jenny  et 
«L  (T|. 

The  final  test  of  the  model,  a  comparison  between  experimental  results  for  stiffness  and  damping  coefficients 
and  the  predictions  of  this  model,  will  be  carried  out  in  Part  2  of  this  paper. 


NOMENCLATURE 


A  Cross-sectional  area  of  control  volume  (L2)\  illustrated  in  figure  (6) 

B  Height  of  labyrinthseal  strip  (L);  illustratedin  figure  (5) 

C  Direct  damping  coefficient  (Ft/L) 

Cr  Nominal  radial  clearance  (L);  illustrated  in  figure  (5) 

Dk  Hydraulic  diameter  of  cavity  (L);  introduced  in  equation  (6) 

H  Local  radial  clearance  (L) 

K  Direct  stiffness  coefficient  (F/L) 

L  Pitch  of  seal  strips  (L);  illustrated  in  figure  (5) 

NT  Number  of  seal  strips 
NC=NT-1  Number  of  cavities 
P  Pressure  {F/L2) 

R  Gas  constant  (L2/Tt2) 

Rs  Radius  of  control  volume  (L);  illustrated  in  figure  (5) 

Rtu  Surface  velocity  of  rotor  (L/t) 

T  Temperature  (T) 

Tp  Tooth  tip  width  (L);  illustrated  in  figure  (5) 

U  Average  axial  velocity  for  control  volume  (L/t);  in  figure  (2) 

W  Average  circumferential  velocity  for  control  volume  (L/t);  illustrated  in  figure  (2) 

W0  Average  circumferential  velocity  in  the  interface  between  control  volumes  I  and  H  (L/t);  introduced 
equation  (18) 

a,b  Radial  seal  displacement  components  due  to  elliptical  whirl  (L);  introduced  in  equation  (51) 
ar,  as  Dimensionless  length  upon  which  shear  stress  acts;  introduced  in  equation  (3)  and  (4) 
c  Cross  coupled  damping  coefficient  (Ft/L);  in  equation  (57) 
e0  Displacement  of  the  seal  rotor  from  centered  position  (L) 
k  Cross  coupled  stiffness  coefficient  (F/L);  in  equation  (57) 
m  Leakage  mass  flow  rate  per  circumferential  length  (M/Lt) 
mr,  nr,  ms,  ns  Coefficients  for  friction  factor;  in  equation  (3) 
t  Time  (t) 

v  Total  velocity  (L/t);  introduced  in  equation  (11) 
w  Shaft  angular  velocity  (l/t) 
p  Density  of  fluid  (M/L3) 
v  Kinematic  viscosity  ( L2/t ) 

<  =  t-o/Cr  Eccentricity  ratio 

t  Turbulent  viscosity  (Ft/L2);  introduced  in  equation  (19) 

7  Ratio  of  specific  heats 

SubteripU 

o  Zeroth-order  component 

1  First-order  component,  control  volume  I  value 

2  Control  volume  II  value 
i  i-th  chamber  value 

j  Value  along  the  dividing  streamline 
x  X -direct  ion 
y  Y-direction 
r  Reservoir  value 
s  Sump  value 
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APPENDIX  A:  GOVERNING  EQUATIONS  FOR  TEETH-ON-STATOR  SEAL 
Reduced  Equation* 

The  main  difference  between  the  teeth-on-stator  equation*  and  the  teeth-on-rotor  equation*  occur*  in  the 
momentum  equation*.*  The  ahear  atreaae*  acting  on  control  volume  I  are  now  the  rotor  (hear  streu,  rr,  and 
the  free  shear  atreaa,  ty.  Similarly,  the  shear  stresses  acting  on  control  volume  II  are  now  the  stator  shear 
stress,  r„  and  the  free  shear  stress,  ry.  These  differences  are  evident  in  the  reduced  form  of  the  continuity 
and  circumferential  momentum  equations  given  below: 

Continuity  I 

dpAi  .  dpW,A,  ,  ^  8pA2  dpW2A2  _n  , 

—  +  ^7i5r+mi+i_m'  +  “ar  +  _R^r-0  {Al) 


Momentum  I 


Momentum  II 


aw,  Pw,A,aw,  (dpA 3  aPw2A2\ 
pAl^  +  ~E^'ir  +  \~aT  +  ~R^dr) 

A  BP 

+  -  W'l.-l)  =  +  T)iLi  -  friariLi 

.  dw2  .  pW2A2avv3  .  ,  a^3A2\,r„  „ 

pA7~aT  +  "ST  +  \~aT  +  (lv* " 

A2  dPi 

~  ~R*2  86  TiiLi  ,ia>iL' 


where  as  and  ar  are  defined  as: 


as,  =  (2B,-  ■+■  Li)/ Li  or<  =  1.0 


The  rotor  shear  stress  in  the  circumferential  direction  is  now  defined  using  the  smaller  hydraulic  diameter 
and  the  velocity  components  of  control  volume  I. 

ty  «  ^py/(R;w-W,y-i-Uf(R,lU  -W,)nr  *SEE (A4) 

where  Dfcl<  is  the  hydraulic  diameter  of  C.V.  I,  defined  by: 

Dun  =  2  CnLiKCnLi)  (A5) 

Similarly,  the  stator  shear  stress  in  the  circumferential  direction  is  now  defined  using  the  larger  hydraulic 
diameter  and  the  velocity  components  of  control  volume  II. 


r,  =  \p<Jwl  +  UZW2n*  ^ 


s/Wl  +  U?Dh 


where  £>j,3<  is  the  hydraulic  diameter  of  C.V.  II,  defined  by: 

Dk2i  =  2  BiLiKBiLi)  (A7) 

The  definition  of  the  free  shear  stress  remains  the  same.  However,  since  CFD  results  were  only  available 
for  the  teeth-on-rotor  configuration,  the  sensitive  mixing  length  ratio,  t/b,  may  change  for  a  teeth-on-stator 
seal. 


Zeroth-Order  Equations 
Continuity 


l  = 


Momentum  I  and  II 


™o(Wioi  —  Wioi-i)  —  (fyoi  —  (i49) 

T} oiU  =  TroiauLi  (>410) 

Firtt-Order  Equations 

The  first-order  equations  remain  exactly  the  aame  aa  before.  Since  change*  were  made  in  the  location*  and 
definition*  of  the  rotor  and  atator  shear  stress  terms,  the  following  changes  in  the  coefficients  of  the  first-order 
equations  are  necessary: 
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APPENDIX  B:  DEFINITION  OF  FIRST-ORDER  COEFFICIENTS 
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APPENDIX  C.  SEPARATION  OF  GOVERNING  EQUATIONS 

Continuity: 
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DEFINITION  OF  MATRIX  ELEMENTS 


Ai-i  MATRIX 

Ai,7  =  Aa.i  =  As,4  =  A«,s  =  Gg 
A®3  =  A«,i  =  A  74  =  A*. 3  =  Xt 
A®,®  =  A*,  ®  =  At,!  =  A*,7  =  ~ "»o 
Ao.a  =  Aio.i  =  Aj  jt4  =  An,®  =  Yt 

The  remaining  element!  are  aero. 

A,  MATRIX 
Ai.i  =  — Aa,2  =  Giu  -f  G2 
Aj.j  =  —  Ae.4  =  Giu  +  Ga 

A  1,2  =  A  2,1  =  Aj,4  =  A4, 3  =  G® 

A®, 2  =  A®,1  =  At,4  *S  A*, 3  =  Xi 

A®,i  =  — A«,a  =  X*u  +  Xa+  - 

,  A  Y  •  y  •  ^Va, 
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A®,®  =  A®,®  =  At,i  —  As,t  =  X® 

Ap.e  =  Aio,s  =  An,*  =  Aia,7  =  Yi 
Ai,p  =  As, ii  =  —  Aa.io  —  -A*, 12  =  G® 

A  A  —  A  —  A  —  XjPcx 
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Ag.io  =  Aio.p  =  An.ia  =  Aia,n  =  V® 

The  remaining  elements  are  sero. 

A,+i  MATRIX 
A  1,2  =  Aa.i  —  A3, 4  =  A4,3  =  G  7 
The  remaining  elements  are  sero. 
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APPENDIX  D:  THE  THEORY  OF  JENNY  ET  AL.  (7} 

Th«  theory  of  Jenny  et  aL  (7j  hu  shown  consistently  good  agreement  with  measured  test  results  [25]  in 
predictions  of  croea-coupled  stiffness  and  direct  damping.  The  author  had  hoped  to  program  their  solution 
and  make  direct  comparison  to  the  present  theory;  however,  as  outlined  below,  unresolvable  difficulties  arose 
in  deriving  the  published  equations  of  (7j. 

The  theory  of  Jenny  et  al.  [7]  was  derived  for  the  *box-in-a-box*  control  volume  configuration  illustrated 
in  figure  5.  Thus,  a  direct  comparison  of  their  equations  with  those  presented  in  this  report  is  not  feasible. 
However,  a  review  of  the  development  of  their  governing  equations  is  of  interest. 

The  following  convention  will  be  used  for  the  control  volumes  in  figure  5:  the  large  control  volume  is  control 
volume  I  and  the  small  control  volume  is  control  volume  II.  The  continuity  equations  for  the  control  volumes 
shown  in  figure  5  are: 

Continuity  1: 


dpW^A? 

RedS 


dpWiAi  dp[Ax  +  A7) 


+  m,+1  -  m,  =  0 


Continuity  II: 


dpV/^Ai  dpAj 


The  following  assumptions  are  used  by  Jenny  et  al.  [7]  to  simplify  equations  (Jl)  and  (J2): 

a)  the  flow  is  incompressible  ( p  =  constant), 

b)  riv+i  =  rhi,  and 

c)  the  area  of  the  control  volume  II  is  constant. 

The  first  assumption  seems  questionable,  since  this  is  a  compressible  flow  solution  and  quite  often  the  flow 
in  a  labyrinth  seal  achieves  Mach  1  at  the  exit.  Assumption  (b)  is  a  valid  assumption  for  the  seroth-order, 
steady  flow  solution,  but  it  is  questionable  for  the  first-order,  unsteady  flow  solution  for  an  orbiting  rotor. 
Using  the  chain  rule  for  the  expansion  of  partial  derivatives  and  the  above  assumptions,  equations  (Jl)  and 
(J2)  become: 

Continuity  I: 


Continuity  II: 


dW-2 

PM-Zi - R»™ri  -  0 

ov 


The  equations  given  by  Jenny  et  al.  [7]  are: 
Continuity  I: 


aw?  aw,  acr  8(a1  +  a3) 

+  Al~df-  WlL^f  9 at  ~ 0 


Continuity  II: 


,  aw7  . 

pAi~df  ~  rn"  =  0  l-76) 

The  difference  between  equations  (J3)  and  (J5)  is  in  the  sign  of  the  third  and  fourth  terms.  The  second  and 
third  terms  in  equations  (J3)  and  (J5)  originate  from  the  same  partial  derivative,  but  have  opposite  signs. 
The  author  could  not  arrive  at  the  same  conclusion  using  the  chain  rule.  The  difference  between  equations 
(J4)  and  (J6)  is  the  radius,  Rs,  in  the  second  term.  This  may  or  may  not  be  a  problem  since  the  radial  mass 
flow  term,  m„,  is  not  defined  by  Jenny  et  al.  |7j. 


The  author  Agreed  with  the  derivAtion  of  the  momentum  equations  for  the  control  volumes  shown  in  figure 
S  except  for  the  Aforementioned  Assumptions  And  the  following  discrepancies: 

(a)  the  axial  velocity  component  is  incorporated  in  the  definition  of  the  stator  wall  shear  stress,  but 
neglected  in  the  definition  of  the  Reynold's  number  which  is  used  to  calculate  the  friction  factor  term  in  the 
shear  stress  relation. 

(b)  the  perturbation  of  the  friction  factor  is  ignored.  This  term  has  been  shown  (24]  to  be  important  in 
the  solution  for  rotordynamic  coefficients. 

(c)  the  leakage  equation  is  a  global  leakage  equation.  This  means  that  local  perturbations  pertaining  to 
a  cavity  can  not  be  found  from  this  equation.  Jenny  et  al.  (7]  perturb  this  global  equation  for  clearance. 

(d)  the  carryover  coefficient  definition  used  in  the  leakage  equation  is  a  global  equation  and  cannot  be 
perturbed. 

(e)  the  flow  coefficient  used  in  the  leakage  equation  was  obtained  from  a  plot  of  empirical  data.  No 
explanation  was  given  for  the  method  used  to  obtain  the  dervatives  of  the  flow  coefficient  used  in  the 
perturbation  equations. 

The  aforementioned  problems  prevented  the  author  from  obtaining  a  solution  based  on  the  theory  of  Jenny 
et  al.  (7j.  Regrettably,  no  direct  comparison  between  it  and  the  theory  presented  in  this  paper  was  possible. 
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Fig.  2  Control  volumes  with 
geometric  boundary. 


Fig.  4  Control  volumes  with 
dividing  streamline 
boundary . 
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Fig.  3  "Box-in-a-box"  model 
of  Jenny  et  al.  7  . 
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Fig.  5  A  typical  labyrinth 
cavity . 


THEORY  VERSUS  EXPERIMENT  FOR  THE 
ROTORDYNAMIC  COEFFICIENTS  OF  LABYRINTH  GAS  SEALS: 
PART  H  -  A  COMPARISON  TO  EXPERIMENT1 
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SUMMARY 

An  experimental  teat  facility  ie  need  to  meaaure  the  leakage  and  rotordynamic  coefficients  of  teeth-on- 
rotor  and  teeth-on-stator  labyrinth  gas  seals.  The  teat  results  are  presented  along  with  the  theoretically 
predicted  values  for  the  two  seal  configurations  at  three  different  radial  clearances  and  shaft  speeds  to  16,000 
cpm.  The  test  results  show  that  the  theory  accurately  predicts  the  cross-coupled  stiffness  for  both  seal 
configurations  and  shows  improvement  in  the  prediction  of  the  direct  damping  for  the  teeth-on-rotor  seaL 
The  theory  fails  to  predict  a  decrease  in  the  direct  damping  coefficient  for  an  increase  in  the  radial  clearance 
for  the  teeth-on-stator  seal 

INTRODUCTION 

Part  1  of  this  paper  presented  a  new  two-control- volume  analysis  to  predict  the  rotordynamic  coefficients  for 
labyrinth  gas  seals.  This  part  (Part  2)  of  the  paper  provides  a  comparison  of  the  predictions  of  the  analysis 
from  Part  1  to  the  new  teat  results  for  six  “see-through*  labyrinth  gas  seals,  as  shown  in  figure  1,  three  with 
teeth  on  the  rotor  and  three  with  teeth  on  the  stator.  The  design,  development,  and  operation  of  the  test 
apparatus  and  facility,  which  have  been  developed  to  measure  the  leakage  and  rotordynamic  coefficients  of 
annular  gas  seals,  has  been  described  by  Childs  et  al.  |l].  The  apparatus  described  in  [l]  was  limited  to 
a  top  shaft  speed  of  8,000  cpm.  The  apparatus  has  since  been  redesigned  to  operate  at  shaft  speeds  up  to 
16,000  cpm.  A  complete  discussion  of  the  redesign  of  the  apparatus  can  be  found  in  Elrod  and  Childs  [2]. 

As  described  in  (l],  the  rotordynamic  coefficients  for  a  gas  seal  are  defined  by  the  following  linearised  force- 
displacement  model. 


(1) 


Where  (X,Y)  define  the  motion  of  the  seal’s  rotor  relative  to  its  stator,  (F»,  F9)  are  the  components  of  the 
reaction  force  acting  on  the  rotor,  and  (K„,  “d  (Cmm,Cn,Cmv,Ctz)  are  the  stiffness  and 

damping  coefficients  respectively.  Equation  (1)  applies  for  small  motion  of  the  rotor  about  an  arbitrary 
eccentric  position.  For  small  motion  about  a  centered  position,  the  following  simpler  mode]  applies. 


{/>}-[-*  £]{y}+[-c  c]{£}  (2) 

Although  the  test  apparatus  has  the  capability  of  separately  identifying  the  eccentric-position  rotordynamic 
coefficients  of  equation  (l),  the  results  presented  here  are  for  the  centered-position  case  on]y. 

PREVIOUS  EXPERIMENTAL  PROGRAMS 

A  limited  amount  of  experimental  data  has  been  published  to  date  on  the  determination  of  stiffness  and 
damping  coefficients  for  labyrinth  gas  seals.  The  first  published  results  for  stiffness  coefficients  were  those 
of  Wachter  and  Benckert  [3,4,5].  They  investigated  the  following  three  types  of  seals:  a)  teeth-on-stator, 
b)  teeth  on  the  rotor  and  stator,  and  c)  teeth  on  the  stator  and  steps  or  grooves  on  the  rotor.  These 
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results  were  limited  in  that  the  pressure  drop  was  small,  much  of  the  data  were  for  nonrotating  seals,  no 
data  were  presented  for  seals  with  teeth  on  the  rotor,  the  rotor  speed  was  limited,  and  tests  where  rotation 
and  inlet  tangential  velocity  existed  simultaneously  were  very  scarce.  The  next  investigation  was  carried 
out  by  Wright  [8],  who  measured  an  'equivalent*  radial  and  tangential  stiffness  for  single-cavity  seals  with 
teeth  on  the  stator.  Although  for  a  very  limited  and  special  case,  Wright's  results  do  give  insight  into  the 
effect  of  pressure  drop,  convergence  or  divergence  of  the  clearance,  and  forward  or  backward  whirl  of  a 
seal.  These  results  could  be  reduced  to  direct  and  cross-coupled  stiffness  and  damping,  hence,  they  are  the 
first  published  damping  coefficients  for  teeth-on-stator  labyrinth  seals.  Brown  and  Leong  [7,8]  investigated 
the  same  seal  configurations  as  Wachter  and  Benckert,  in  an  effort  to  verify  and  extend  their  work.  Their 
results  include  variations  of  pressure,  geometry,  rotor  speed,  and  inlet  tangential  velocity.  Although  the 
investigation  was  extensive,  the  published  results  are  limited  because  of  the  lack  of  information  concerning 
operating  conditions  for  the  various  tests.  Childs  and  Scharrer  [0]  investigated  geometrically  similar  teeth- 
on-rotor  and  teeth-on-stator  labyrinth  gas  seals  for  stiffness  and  damping  coefficients  up  to  speeds  of  8000 
cpm.  Kanemitsu  and  Ohsawa  [10]  investigated  multistage  teeth-on-stator  and  interlocking  labyrinth  seals 
up  to  speeds  of  2400  cpm.  They  measured  an  effective  radial  and  tangential  stiffness  while  varying  the  whirl 
frequency  of  the  rotor.  These  data  could  be  reduced  to  stiffness  and  damping  coefficients.  Hisa  et  ai  (11} 
investigated  teeth-on-stator  seals  with  2-4  teeth  and  a  teeth-on-stator  seal  with  steps  on  the  rotor  up  to 
speeds  of  8000  cpm.  These  data  only  included  static  tests  for  direct  and  cross-coupled  stiffness  using  steam. 

The  most  extensive  comparison  of  analytical  predictions  and  experimental  results  was  carried  out  by  Scharrer 
[12]  using  the  theory  of  Childs  and  ScharTer  [13]  and  the  results  of  Childs  and  Scharrer  [9].  This  comparison 
showed  that  the  theory  [13]  predicts  cross-coupled  stiffness  reasonably  well,  but  underpredicts  direct  stiffness, 
direct  damping  and  cross-coupled  damping. 

In  reviewing  previous  experimental  programs,  there  is  a  dear  need  for  a)  more  extensive  testing  of  teeth- 
on-rotor  labyrinth  seals  and  b)  experimental  results  showing  the  effect  of  radial  clearance  change  on  direct 
damping  coefficients.  This  paper  addresses  these  points  in  addition  to  evaluating  the  new  analysis  presented 
in  Part  1  of  this  paper  by  comparison  to  the  new  test  results. 

TEST  APPARATUS  AND  FACILITY 
Introduction 

The  test  results  reported  here  were  developed  as  a  part  of  an  extended,  joint  NASA-USAF  funded  research 
program  for  annular  gas  seal  studies.  Tests  were  carried  out  on  six  *  see-through’  labyrinth  seals,  three 
with  teeth  on  the  rotor  and  three  with  teeth  on  the  stator,  each  with  different  radial  clearances.  The  test 
program  had  the  objective  of  examining  the  effects  of  a  change  in  radial  clearance  on  the  leakage  and  stability 
performances  of  a  teeth-on-stator  and  a  teeth-on-rotor  labyrinth  seal.  Air  is  the  test  fluid. 

Test  Apparatus 

The  rotor  shaft  is  suspended  pendulum-fashion  from  an  upper,  rigidly  mounted  pivot  shaft,  as  shown  in  figure 
2.  This  arrangement  allows  a  side-to-side  (horisontal)  motion  of  the  rotor.  A  cam  within  the  pivot  shaft 
allows  vertical  positioning  of  the  rotor.  The  rotor  is  both  positioned  and  excited  horisontally  by  a  hydraulic 
shaker  head  which  acts  on  the  rotor-shaft  bearing  housing  and  works  against  a  return  spring  mounted  on  the 
opposite  side  of  the  bearing  housing.  The  design  of  the  test  rig  permits  the  installation  of  various  rotor /stator 
combinations.  The  stator  is  supported  in  the  test  section  housing  by  three  pieso-  electric  quarts  load  cells 
in  a  trihedral  configuration.  The  test  apparatus  measures  the  reaction-force  components  and  relative  seal 
stator  motion.  The  harmonic  components  of  the  motion  and  force  signals  are  used  to  identify  the  stiffness 
and  damping  coefficients.  Different  seal  stator  designs  are  obtained  by  the  use  of  inserts. 

The  dimensions  and  pertinent  data  for  each  seal  configuration  are  given  in  table  1.  For  the  remainder  of 
this  paper,  the  seals  will  be  referred  to  as  seal  1,  seal  2  and  seal  3,  as  given  in  table  1,  in  addition  to  their 
respective  configuration.  The  smooth  and  labyrinth  inserts  used  for  the  0.4mm  (O.OlOin.)  clearance  seals 
are  shown  in  figure  3.  The  labyrinth  tooth  detail  for  both  rotor  and  stator  is  shown  in  figure  4. 


Table  1.  Dimension*  of  seals  tested  in  this  study 
Teeth-on-rotor 


Teeth-on-stator 


Seal  1 

Radius,  mm 

72.5 

75.6 

Length,  mm 

50.8 

50.8 

Clearance,  mm 

0.3 

0.33 

Number  of  teeth 
Seal  2 

16 

16 

Radius,  mm 

72.5 

75.6 

Length,  mm 

50.8 

50.8 

Clearance,  mm 

0.4 

0.4 

Number  of  teeth 
Seals 

16 

16 

Radius,  mm 

72.5 

75.6 

Length,  mm 

50.8 

50.8 

Clearance,  mm 

0.55 

0.6 

Number  of  teeth 

16 

16 

Teet  Variable $ 

When  shaking  about  the  centered  position,  the  Dynamic-Seal-Apparatus  is  capable  of  controlling  the  fol¬ 
lowing  three  independent  variables: j>r***sr*  ratio,  rotor  speed  and  inlet  circumferential  ve/octty.Two  shake 
frequencies,  56.8  Hs  and  74.6  Hi,  were  used  during  testing  with  essentially  the  same  results.  The  results  to 
be  presented  were  obtained  by  shaking  at  74.6  Hs  at  an  amplitude  between  0.076  mm  and  0.1  mm.  The 
actual  test  points  for  each  of  these  three  independent  variables  are  shown  in  table  2.  When  reviewing  the 
following  figures,  table  2  should  be  consulted  for  the  definitions  of  all  symbols  used. 

Table  2.  Definition  of  symbols  used  in  figures 

Supply  pressure 

Rotor  speeds 

Inlet  circumferential  velocities 

1  3.08  bar 

1  3000  cpm 

1  High  vel.  against  rotation 

2  4.46  bar 

2  6000  cpm 

2  Low  veL  against  rotation 

3  5.84  bar 

3  9500  cpm 

3  Zero  circumferential  vel. 

4  7.22  bar 

4  13000  cpm 

4  Low  vel.  with  rotation 

5  8.22  bar 

5  16000  cpm 

5  High  vel.  with  rotation 

The  reservoir  pressures,  as  measured  upstream  of  the  flowmeter,  are  given  in  table  2.  These  values  differ 

from  the  actual  inlet  pressure 

because  of  frictional  losses  and  an 

acceleration  of  the  fluid  due  to  the  inlet 

guide  vanes.  No  tests  could  be  run  at  sero  pressure  difference,  since  a  small  pressure  difference  is  necessary 

to  keep  the  rotor  from  «hifti«v 

axially  and  rubbing  the  inlet  guide 

vanes.  Similarly,  no  sero  rotor  speed  tests 

were  run,  since  rotor  rotation  was  necessary  to  prevent  damage  to  the  thrust  bearing  during  shaking. 

The  inlet  circumferential  velocities  are  given  in  figures  5,6  and  7  as  a  function  of  pressure  ratio.  The  teeth- 
on-rotor  results  are  on  the  left  and  the  teeth-on-stator  results  are  on  the  right.  This  convention  will  be  used 
for  all  of  the  results  presented  in  this  paper.  The  figures  show  that  inlet  circumferential  velocity  remains 
fairly  constant  over  the  pressure  ratios  tested.  There  were  five  test  points  for  inlet  circumferential  velocity; 
two  positive,  two  negative,  and  one  at  sero.  The  sero  inlet  circumferential  velocity  point  corresponds 
to  the  x-axis  in  the  figures  5,6  and  7.  The  negative  numbers  shown  in  the  figures  mean  that  the  inlet 
circumferential  velocity  was  opposed  to  the  direction  of  rotor  rotation.  The  positive  numbers  mean  that  the 
inlet  circumferential  velocity  was  in  the  same  direction  as  rotor  rotation.  The  two  different  magnitudes  of 
inlet  circumferential  velocity,  for  each  direction,  correspond  to  the  different  inlet  guide  vane  geometries,  as 
discussed  in  (lj.  The  ratio  of  inlet  circumferential  velocity  to  rotor  surface  velocity  (inlet  circumferential 
velocity  ratio),  ranged  from  about  -6  to  about  6.  Although  the  larger  numbers  are  practically  unrealistic, 
they  do  give  insight  into  the  effects  of  inlet  circumferential  velocity  that  would  have  otherwise  gone  unnoticed. 


Normalized  Parameter! 

Before  the  teste  described  herein  were  performed,  the  TAMU  gas  seal  test  apparatus  was  modified  as  de¬ 
scribed  by  Elrod  and  Childs  (2]  to  allow  operation  at  running  speeds  up  to  16,000  cpm.  As  expected, 
subsequent  tests  revealed  a  dependence  of  the  rotor  diameter  on  running  speed  due  to  inertia  and  thermal 
effects.  The  rotor  growth  data,  shown  in  table  S,  were  obtained  from  eddy  current  motion  probes  positioned 
at  the  midspan  of  the  seal.  Thus,  as  the  rotor  turns  faster,  the  forces  in  the  seal  are  affected  not  only  by  the 
increased  surface  speed  of  the  rotor  but  also  by  a  change  in  clearance.  Theoretically,  normalisation  would 
collapse  the  data  and  make  the  presentation  simpler  and  more  straight  forward.  This  was  not  the  case  with 
the  labyrinth  seals  tested  in  this  study.  The  failure  of  the  normalisation  is  detailed  by  Scharrer  [14]. 

Table  3.  Growth  of  rotor  with  rotational  speed 

Rotor  speed  (cpm)  Diametrical  growth  (mm) 

3000  0.01 

6000  0.02 

9500  0.03 

13000  0.05 

16000  0.11 

Dynamic  Results 

For  a  circular  orbit  of  amplitude  A,  the  resultant  radial  and  tangential  forces  developed  by  the  seal  model 
of  equation  (2)  are  illustrated  in  figure  8  and  are  defined  by 

-Ft  l A  =  K  +  cu 
Ft) A  =  k  —  Cu 

From  a  stability  standpoint,  the  destabilising  tangential  force,  Ft,  is  of  most  interest.  The  destabilising 
influence  comes  from  the  cross-coupled  stiffness,  k,  and  the  stabilising  influence  comes  from  the  direct 
damping,  C.  The  radial  force  usually  has  little  influence  on  stability,  except  in  rare  cases  involving  multistage 
*  back-to-back*  centrifugal  compressors  with  midapan  seals  where  large  negative  direct  stiffness  values  may 
reduce  the  natural  frequencies.  Since  the  focus  of  this  study  was  on  stability,  the  cross-coupled  stiffness  and 
direct  damping  results,  which  have  the  most  influence,  will  be  presented  first.  The  direct  stiffness  will  follow. 

Relative  Uncertainty 

Before  proceeding  with  the  results,  a  statement  must  be  made  concerning  the  uncertainty  present  in  the 
experimental  results.  Using  the  method  described  by  Holman  [15],  the  uncertainty  in  the  dynamic  coefficients 
can  be  determined.  The  uncertainty  in  the  force,  excitation  frequency,  and  displacement  measurements  are 
0.89  N  (0.2  lb),  0.13  Hs  ,  and  0.0013  mm  (0.05  mils),  respectively.  The  resulting  calculated  uncertainty  is  7 
N/mm  (40  lb/in)  for  the  stiffness  coefficients  and  0.0875  N-s/mm  (0.5  lb-s/in)  for  the  damping  coefficients. 
Since  the  measured  cross-coupled  damping  results  were  rarely  greater  than  the  uncertainty,  test  results  are 
not  provided  here  for  this  parameter. 

TEST  RESULTS 

It  might  seem  obvious,  since  this  paper  evaluates  the  effect  on  seal  performance  of  varying  the  radial  seal 
clearance,  that  the  data  should  be  presented  as  a  function  of  clearance  (clearance  being  the  x-axis).  However, 
since  the  inlet  circumferential  velocity  is  directly  dependent  on  seal  leakage  and  the  seals  leakage  a  different 
rates  due  to  differing  cross-sectional  areas,  the  inlet  circumferential  velocity  test  points  for  seal  1  differs  from 
those  of  seals  2  and  3.  This  is  a  problem  because  the  rotordynamic  coefficients  are  very  sensitive  to  the  inlet 
circumferential  velocity.  Therefore  the  dynamic  data  will  be  presented  as  a  function  of  inlet  circumferential 
velocity  ratio  at  one  pressure  and  one  rotor  speed.  In  the  following  figures,  the  solid  lines  represent  the  test 
results  and  the  broken  lines  represent  the  analytical  predictions. 

Cross-coupled  Stiffness  Comparison 

Figure  9  shows  a  comparison  of  experimental  and  theoretical  results  for  cross-coupled  stiffness  versus  inlet 
circumferential  velocity  ratio  for  the  three  seal  clearances  of  table  1.  The  figure  shows  that  the  theory  does  a 
good  job  of  predicting  the  cross-coupled  stiffness  for  both  the  teeth-on-rotor  and  teeth-on-stator  seals,  at  low 
rotor  speeds.  The  figure  also  shows  that  there  is  no  consistent  trend  in  the  cross-coupled  results  for  a  change 
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in  clearance.  Figure  10  ahowa  a  companion  of  experimental  and  theoretical  remits  for  cross-coupled  stiffnesa 
verms  rotor  speed  for  seal  1  (minimum  clearance  seal)  of  table  1  at  the  inlet  pressures  of  table  2.  The  figure 
shows  that  the  theory  correctly  predicts  a  sharp  upturn  in  cross-coupled  stiffnesa  at  the  higher  rotor  speeds, 
for  the  teeth-on-rotor  seal.  This  speed  sensitivity  was  not  evident  in  the  low  speed  results  of  [Oj  nor  in  the 
remits  for  seals  2  and  $.  The  figure  also  shows  that  the  cross-coupled  stiffnesa  for  the  teeth-on-stator  seal 
decreases  as  rotor  speed  increases.  This  decrease  in  cross-coupled  stiffness  for  an  increase  in  rotor  speed  was 
also  reported  by  Hisa  et  aL  (11]  for  their  teeth-on-stator  seals.  The  theory  shows  the  same  decrease  until  the 
higher  rotor  speeds  are  reached,  then  h  shows  a  sharp  upturn.  This  upturn  may  appear  in  the  experimental 
data  at  higher  rotor  speeds.  The  same  decrease  in  cross-coupled  stiffness  with  an  increase  in  rotor  speed  was 
evident  in  the  results  for  seals  2  and  S. 

Direct  Damping  Comparison 

Figaro  11  shows  a  comparison  of  experimental  and  theoretical  direct  damping  versus  inlet  circumferential 
velocity  ratio  at  the  three  clearances  of  table  1.  The  figure  shows  that  the  direct  damping  coefficient  for  a 
teeth-on-rotor  seal  increases  as  clearance  increases.  The  theory  correctly  predicts  this  trend  but  underpredicts 
the  magnitude  of  the  coefficient  by  SO  percent.  The  figure  also  shows  that  the  direct  damping  coefficient 
for  a  teeth-on-stator  seal  decreases  as  clearance  increases.  The  theory  shows  the  opposite  trend.  Figure  12 
shows  a  comparison  of  experimental  and  theoretical  direct  damping  versus  rotor  speed  for  seal  1  of  table 
1  at  the  inlet  pressures  of  table  2.  The  figure  shows  that  the  test  remits  show  little  or  no  sensitivity  to 
rotor  speed,  for  either  seal  configuration,  while  the  theory  shows  a  sharp  upturn  at  the  higher  speeds.  This 
upturn  may  be  evident  in  future  test  results  at  higher  speeds.  The  results  for  seals  2  and  3  showed  the  same 
insensitivity  to  rotor  speed. 

Direct  Stiffness  Comparison  — 

Figure  13  shows  a  comparison  of  experimental  and  theoretical  direct  stiffness  versus  inlet  circumferential 
velocity  for  the  seal  clearances  of  table  1.  The  figure  shows  that,  for  both  seal  configurations,  the  direct 
stiffness  coefficient  is  negative  and  increases  as  radial  seal  clearance  increases.  One  would  expect  sero 
direct  stiffnesa  values  at  sufficiently  large  clearances.  The  theory  predicts  this  same  trend,  for  both  seal 
configurations.  Figure  13  also  shows  that  the  theory  underpredicts  the  direct  stiffness  magnitudes  at  high 
rotor  speeds.  Figure  14  show  a  comparison  of  experimental  and  theoretical  direct  stiffness  versus  rotor  speed 
for  seal  1  of  table  1  at  the  inlet  pressures  of  table  2.  The  figure  shows  that  the  test  results,  for  both  seal 
configurations,  show  little  or  no  sensitivity  to  rotor  speed.  The  figure  also  shows  that  the  theory  is  overly 
sensitive  to  rotor  speed.  The  insensitivity  to  rotor  speed  was  also  evident  in  the  test  results  for  seals  2  and 
3. 

Stability  Analysis 

One  of  the  main  objectives  of  this  test  program  was  an  evaluation  of  the  effect  on  seal  performance  of  varying 
the  radial  seal  clearance.  A  comparison  of  the  stability  of  the  two  seal  configurations  at  the  different  radial 
clearances  satisfies  that  objective.  Since  a  direct  comparison  of  the  coefficients  of  the  two  seals  does  not 
show  any  clear  stability  advantage  or  the  overall  effect  of  a  change  in  the  radial  seal  clearance,  another 
method  of  comparison  must  be  used.  One  method  in  which  the  dynamic  coefficients  of  the  two  seals  can  be 
directly  compared  is  through  their  respective  non-dimensional  whirl  frequency  ratios.  Whirl  frequency  ratio 
is  defined  by 


Whirl  frequency  ratio  =  k/Cw  =  0/ 

where  u  is  the  running  speed,  and  0/  is  the  ratio  of  the  destabilising  influence  of  the  cross-coupled  stiffness 
and  the  stabilising  influence  of  direct  damping.  fVom  a  stability  viewpoint,  a  minimum  whirl  ratio  is 
desirable.  Figure  15  shows  the  whirl  frequency  ratio  versus  inlet  circumferential  velocity  ratio  for  the  seal 
clearances  of  table  1.  For  teeth-on-rotor  seals,  the  figure  shows  that,  as  the  clearance  increases,  the  seal 
becomes  more  stable.  For  teeth-on-stator  seals  the  opposite  is  true;  as  clearance  increases  the  seal  becomes 
leas  stable,  for  the  positive  inlet  circumferential  velocity  case.  The  figure  also  shows  that  the  teeth-on-stator 
seals  are  more  stable  that  the  teeth-on-rotor  seals  for  the  positive  inlet  circumferential  velocity  ratio,  as  was 
found  previously  [9], 
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Companion  to  tkeory  of  / 13} 

Figure*  10  ud  17  provide  a  brief  comparison  of  the  present  theory  to  the  theory  of  Childs  and  Schsrrer  (IS], 
Figure  10  shows  cross-coupled  stiffness  venns  pressure  ratio  for  a  teeth-on-rotor  labyrinth  seal  at  10,000 
cpm,  and  demonstrates  that  the  present  theory  follows  the  experimental  data  mare  closely  than  the  former 
theory  {IS].  Similarly,  .figure  17  shows  that  the  present  theory  also  follows  the  experimental  data  for  direct 
damping  more  closely  than  the  former  theory  (IS]. 

CONCLUSIONS 

Test  results  have  been  presented  for  stiffness  and  damping  coefficients  for  six  See-through”  labyrinth  seals, 
three  with  teeth  on  the  rotor  and  three  with  teeth  on  the  stator.  The  seals  were  tested  under  identical 
operating  conditions  to  investigate  the  influence  of  rotor  speed  and  the  effect  of  varying  the  radial  seal 
clearance  on  the  rotordynamic  coefficients.  These  experimental  results  were  compared  to  the  predictions 
from  the  new  analysis  presented  in  Part  1  of  this  paper. 

The  experimental  results  of  the  previous  section  support  the  following  conclusions: 

(1)  For  teeth-oo -rotor  seals,  the  direct  damping  increases  as  clearance  increases;  for  teeth-on-etator 
seals,  the  direct  damping  decreases  as  clear  an  css  increases. 

(2)  Direct  stiffness  is  negative  and  increases  as  clearances  increases,  for  both  seal  configurations.  Cross- 
coupled  stiffness  showed  no  consistent  trend  with  respect  to  clearance  changes. 

(S)  Direct  stiffness  and  direct  damping  show  little  or  no  sensitivity  to  rotor  speed  up  to  16,000  cpm. 
Croes-coupled  stiffness  shows  a  sharp  upswing  at  higher  rotor  speeds,  for  a  teeth-on-rotor  seal.  Cross-coupled 
stiffness  decreases  as  rotor  speed  increases,  for  a  teeth-on-etator  seal 

(4)  As  clearance  decreases,  teeth-on-rotor  seals  become  less  stable  and  teeth-on-stator  seals  become 
more  stable,  for  positive  inlet  circumferential  velocity. 

The  theoretical  results  of  the  previous  section  support  the  following  conclusions: 

(1)  The  theory  correctly  predicts  that  direct  stiffness  is  negative  and  increases  as  clearance  increases,  for 
both  seal  configurations.  The  theory  incorrectly  predicts  an  approximately  quadratic  increase  in  the  direct 
stiffness  magnitude  (becoming  more  negative)  as  speed  increases.  Test  results  show  scant  sensitivity. 

(2)  The  theory  accurately  predicts  an  increase  in  cross-coupled  stiffness  at  high  speeds,  for  a  teeth-on- 
rotor  seaL 

(3)  For  teeth-on-rotor  seals,  the  theory  correctly  predicts  an  increase  in  direct  damping  for  an  increase 
in  clearance.  However,  the  theory  incorrectly  predicts  the  same  trend  for  a  teeth-on-stator  seal. 

(4)  The  theory  incorrectly  predicts  an  approximately  quadratic  increase  in  direct  damping  with  running 
speed.  Test  results  show  no  systematic  change  in  direct  damping  with  running  speed. 

(5)  A  comparison  with  test  results  for  a  teeth-on-rotor  seal  shows  that  the  theory  presented  in  Part  1 
of  this  paper  does  a  better  job  of  predicting  cross-coupled  stiffness  and  direct  damping  than  does  the  theory 
of  Childs  and  Scharrer  (13]. 
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NOMENCLATURE 


A  Seal  orbit  radius  (L);  illustrated  in  figure  8. 

B  Tooth  height  (L);  illustrated  in  figure  1. 

C,c  Direct  and  crow  coupled  damping  coefficient*  (FT/L) 
Cr  Redial  clearance  (L);  illustrated  in  figure  1. 

K,k  Direct  and  cross-coupled  stiffness  coefficients  (F/L) 

F  Seal  reaction-force  (F) 

L  Tooth  pitch  (L);  illustrated  in  figure  1. 

Pr  Seal  inlet  pressure  ( F/L a) 

Rs  Seal  radius  (L);  illustrated  in  figure  1. 

X,Y  Rotor  to  stator  relative  displacement  components  (L) 
0/  Whirl  frequency  ratio 
u  Shaft  angular  velocity  (l/T) 

Subscripts 

i  Value  in  i-th  cavity 
r  Radial  component 
t  Tangential  component 
zjr  Rectangular  coordinate  directions 
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Fig.  10  A  comparison  of  experimental  and  theoretical  cross-coupled 
stiffness  versus  rotor  speed  for  seal  1  and  Inlet 
circumferential  velocity  5.  See  table  2  for  definitions. 
Teeth-on-rotor  (left),  teeth-on-stator  (right). 
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Fig.  14  A  coaparlflon  of  experimental  and  theoretical  direct  stiffness 
versus  rotor  speed  for  seal  1  and  inlet  clrcianferentlal 
velocity  5.  See  table  2  for  definitions.  Teeth-on-rotor 
(left)!  teeth-on-stator  (right). 
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Fig.  15  Whirl  frequency  ratio  versus  inlet  circumferential  velocity 
ratio  at  an  Inlet  pressure  of  3.08  bar  and  rotor  speed  of 
16000  cpra .  Shake  frequency  7*. 6  Hz.  See  table  1  for  seal 
clearance  definitions.  Teeth-on-rotor  (left),  teeth-on-stator 
(right). 
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16  A  comparison  of  experimental  and  theoretical 
results  of  this  report  with  those  of  [13] 
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Fig.  17  A  comparison  of  experimental  and  theoretical 
results  of  this  report  with  those  of  [13  ] 
for  direct  damping. 
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ABSTRACT 


A  numerical  method  employing  a  finite  difference  approach  for  calculating 
the  rotordynamic  force  on  eccentric, whirling,  labyrinth  seals  has  been  developed. 
The  SIMPLER  algorithm  along  with  QUICK  differencing  is  used  to  calculate  the 
flowfield  within  a  seal.  A  modified  bipolar  coordinate  system  accurately  describes 
the  geometry  of  an  eccentric  seal.  The  high  Reynolds  number  k  —  e  turbulence 
model  is  utilized,  which  can  handle  subsonic  compressible  or  incompressible  flows. 
A  three-percent  eccentric  single  labyrinth  cavity  rotating  at  5,000  c.p.m.  was 
investigated  with  three  different  inlet  swirl  conditions,  each  with  and  without  a 
whirl  orbit  frequency  of  2500  c.p.m.  The  fluid  was  air  with  an  inlet  axial  velocity 
near  Mach  0.2.  Detailed  force,  pressure  and  shear  stress  distributions  within  the 
cavity  are  presented.  The  results  indicate  that  the  pressure  component  accounts 
for  99  percent  of  the  rotordynamic  force.  Whirl  seems  to  have  little  effect  on  the 
force,  and  the  downstream  tooth  of  the  cavity  makes  a  very  significant  contribution 
to  this  quantity. 
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radius  of  tooth  periphery 
stator  radius 
clearance 

convection  coefficient 
turbulence  constant 
turbulence  constant 
turbulence  constant 
diffusion  coefficient 

displacement  along  general  orthogonal  coordinate  lines 
rotordynamic  force 
scale  factor 

coordinate  variation  term 
turbulence  kinetic  energy 
turbulence  length  scale 
seal  cavity  pitch 
pressure 

reference  pressure 
pressure  re-defined  by  Eq.  23 
polar  coordinate 
distance 

x 1  velocity  component 

velocity  expressed  in  cartesian  tensor  notation 

velocity  expressed  in  general  orthogonal  coordinate  system 

x2  velocity  component 
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nondimensional  quantity 


Labyrinth  seals  are  usually  used  in  pumps,  compressors,  and  turbines  to  limit 
internal  leakage.  The  labyrinth  seal  was  first  applied  in  a  steam  turbine  near  the 
turn  of  the  century.  The  idea  was  to  create  a  flow  path  between  high  and  low 
pressure  zones  that  would  convert  pressure  head  into  kinetic  energy,  which  was  then 
dissipated.  This  was  accomplished  by  a  series  of  cavities.  The  restriction  between 
each  pair  of  cavities  in  the  flow  path  converts  the  pressure  head  into  mean  flow 
kinetic  energy,  and  the  expansion  in  the  flow  path  dissipates  the  kinetic  energy  by 
viscous  losses.  Figure  1  shows  a  bank  of  generic  labyrinth  cavities  and  the  resulting 
streamlines  due  to  leakage.  Although  seals  are  very  successful  in  limiting  internal 
leakage,  they  are  one  of  the  sources  of  the  self-excited  vibration  of  the  turbomachine 
rotor. 

Self-excited  vibration  is  generally  subsynchronous,  i.e.  the  rotor  whirls  within 
it’s  housing  at  a  frequency  below  the  rotational  frequency.  This  vibration  limits 
operating  speeds  and  can  be  extremely  destructive  to  the  turbomachine.  Although 
the  physics  of  the  shaft  vibration  is  understood,  a  quantitative  analysis  of  the  forces 
which  cause  the  excitation  is  not  complete.  The  cause  of  self  vibration  is  a  net  force 
imbalance  on  the  rotor,  and  this  imbalance  occurs  when  the  rotor  is  displaced  from 
the  geometric  center  of  its  housing. 

Numerous  attempts  to  analyze  the  destabilizing  forces  due  to  leakage  through 
labyrinth  seals  have  been  made.  Analytical,  numerical,  and  experimental  investi¬ 
gations  have  been  conducted.  Most  of  the  analytical  models  are  complex  and  use 
restrictive  assumptions.  Fujikawa  1]  and  Jenny  [2]  used  a  two-control-volume  ana¬ 
lytical  approach,  and  both  used  empirical  data  extensively.  Iwatsubo  [3]  considered 
time  dependent  flow  area  change,  but  neglected  area  change  in  the  circumferential 


direction.  Recently,  Childs  [4,5,6]  improved  Iwatsubo’s  [3]  method  and  obtained 
improved  results.  Kurohashi  [7]  introduced  a  variable  flow  coefficient  for  annular  seals 
and  obtained  reasonable  results  for  the  cross-coupled  force.  The  fluid  mechan¬ 
ics  associated  with  an  eccentric  labyrinth  seal  may  be  too  complex  to  accurately 
model  using  an  analytical  method.  Dietzen  and  Nordman  [8]  used  a  finite  difference 
method  in  conjuction  with  a  perturbation  analysis  to  determine  force  coefficients 
for  annular  seals.  Although  the  method  appears  to  be  promising,  the  author  did 
not  compute  labyrinth  seals.  Baumgartner  [9;  used  analytical  methods  to  simplify 
the  equations,  which  were  then  solved  numerically.  Tam  and  Przekwas  [10]  used 
a  full  three-dimensional  finite  difference  analysis  and  obtained  the  most  detailed 
results  to  date. 

Experimental  measurements  of  the  rotordynamic  coefficients  associated  with 
labyrinth  seals  are  quite  scarce.  Wachter  and  Benckert  [11,12,13)  were  the  first 
to  measure  stiffness  coefficients,  but  much  of  the  results  are  for  nonrotating  seals. 
Wright  [14]  measured  single  cavity  seals  with  teeth  on  the  rotor.  Childs  [4,5,6] 
appears  to  have  obtained  the  most  comprehensive  set  of  measurements.  Rotational 
speed,  inlet  swirl  velocity,  and  pressure  drop  across  the  seals  were  all  varied  in 
these  measurements.  Seals  with  teeth  on  the  rotor  as  well  as  teeth  on  the  stator 
were  used.  Rajakumar  [15]  measured  pressure  in  the  circumferential  direction  and 
measured  forces  at  large  eccentricities. 

COMPUTATIONAL  APPROACH 

For  the  present  study  it  has  been  assumed  that  the  shaft  undergoes  a  circular 
whirl  orbit,  the  center  of  which  coincides  with  that  of  the  housing.  Further,  a 
single  cavity  of  a  multi-cavity  seal  of  teeth-on-rotor  design  as  shown  in  Figure  1  is 
considered.  Note  that  the  periphery  of  the  teeth  and  the  base  of  the  seal  cavity  are 


9 


concentric.  Thus,  cylindrical  coordinates  are  well  suited  to  describe  the  geometry 
within  the  seal  cavity.  However,  in  the  tooth  clearance  region  outside  of  the  cavity, 
bipolar  coordinates  are  used  because  of  the  eccentricity  of  the  shaft.  Specifically, 
it  is  the  modified  bipolar  coordinates  which  were  choosen  because  they  reduce  to 
cylindrical  coordinates  in  the  limit  as  eccentricity  approaches  zero.  Figure  2  shows 
such  a  finite  difference  grid  within  the  gap  region. 

Wood  [16]  used  modified  bipolar  coordinates  to  solve  the  fluid  dynamics 
associated  with  eccentric  rotating  cylinders  in  general.  Consider  two  circular 
cylinders  of  radii  a  and  b.  Let  the  centers  of  the  two  cylinders  be  a  distance  ae 
apart.  If  the  eccentricity  £  and  the  clearance  c  are  defined  as 


ae 

c 

(1) 

and 

c  =  (b  —  a) 

(2) 

then  modified  bipolar  coordinates 

are  defined  by  the  transformation 

+  7) 

1  +  it 

(3) 

z  =  re 16 

(4) 

Cl 

II 

w 

(5) 

where 

7  =  ~2e{(b/a)2  -  1  - 

-  e2  +  y/{b2/a2  -  1  -  e2 )2  -  4E2]'1 

(6) 

and 

,  _  b/a  +  f  -  h 

1  -  (6/flh  -  n 


(7) 


Figure  3  shows  the  geometry  and  coordinate  system  associated  with  modified 
bipolar  coordinates.  Curves  of  constant  p  are  circles  in  the  modified  bipolar 
coordinate  system.  The  inner  circle  is  p  —  1,  and  the  outer  is  p  —  (3.  The  conformal 
transformation  shown  above  provides  a  method  to  transform  any  point  ( p,<f> )  in  the 
modified  bipolar  coordinate  system  into  an  equivalent  point  (r,  0)  in  the  cylindrical 
coordinate  system. 

The  transport  equations  to  be  solved  are  presented  in  Cartesian  tensor  notation. 
Neglecting  body  forces,  the  steady-state  continuity  and  momentum  equations  are 

d(pUt)  _  n 


djpUjUj)  _  dp  dTij 

dx{  dxj  dxi 

where  p,  Ut,  and  p  are  the  time  averaged  density,  velocity  and  pressure.  The 
Reynolds  stresses  are  determined  by  the  k  —  e  turbulence  model.  The  turbulent 
energy  and  energy  dissipation  transport  equations  are 

a 

OX,  OXt  (7k  OX  i 

d(pUte)  d  (P'ff  de  .  e 

-^r  =  dTi{—8i:)  +  -ktc''F-c’2f,i}  (n) 

where  the  stress  is  calculated  using 

,dl h  dU3x  2  ,  3U,C 

T'J  ~  Peff(  dx3  +  dXl)+  3{pk  +  Peff  dxt )6l}  (1“) 

The  production  of  kinetic  energy  and  the  effective  viscosity  are  given  by 

p  =  -(|^)T„  (13) 


Peff  = 


C\pk: 


y_vv 
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The  constants  CM,  Cf],  Ct2,  0k,  and  at  are  given  values  of  .09,1.44,1.92,1.0, 
and  1.217.  This  is  in  accordance  with  the  standard  -  t  model  investigated  by 
Launder  and  Spalding  [17].  For  inlet  Mach  numbers  greater  than  0.3,  the  energy 
equation  must  also  be  employed.  The  stagnation  enthalpy  form 


W) 

dx. 


is  included  where 


d  /T,  d  rrT  „  d 


rH  = 


Static  temperature  is  calculated  using 


„  H-£ 

CD 


The  above  equations  can  be  transformed  into  relations  using  general  orthogonal 
coordinates.  Pope  [18]  suggests  such  a  procedure  using  general  tensor  notation. 
Under  this  procedure  the  equations  retain  their  original  form  and  simplicity  as  much 
as  possible,  and  they  remain  in  a  form  where  each  term  has  physical  significance. 

The  divergence  operator  V(i)  and  coordinate  variation  term  Hx(j)  needed  to 
transform  the  transport  equations  into  general  orthogonal  coordinates  are 


V(i)  = 


hi  d  ]hj 
jh)  dx(i)  h: 


Hr(j)  = 


i  ah, 

h,h;  dxj 


Here  hx  represent  scale  factors  and  |h|  is  the  product  of  the  scale  factors. 
Following  Pope's  18;  transformation  procedure,  the  transport  equations  become 


v(i)>no]  -  o 


Vj»  *V AV  ***-Mm.  IkTJtrJBrLJ**  Ah  'j**  -  “  »  '  J»  **  *  »  *  J.V  M  *  m m*  -*“*  b  "  “  A  ■  *  fc  *“«  »  A  .V  •.  »  .*<  V 


NX 


V(i){pU(i)U(j)  +  T'(ij)}  =  - 


+  Ht{j)[pU(i)V(i)  +  Tm(ii)} 


V(i)\pU(i)iP 


dx(j) 

-H}(i)\pU(i)U(j)  +  T-m 

Peff  dip 


(21) 


=  5, 


(22) 


dx{i)] 

The  last  equation  is  a  general  transport  equation  for  any  quantity  ip.  Note  that 
the  isotropic  component  of  stress  has  been  added  to  the  pressure  giving 


2  2 

P*  =  P  +  ^Pk  +  gMe//V(i)t/(i) 


(23) 


and  r*  has  the  anisotropic  stress  as 

t'(ij)  =  "  V{i)Hi(j)  -  U(j)Hj(i)  +  2 U(l)Ht(l)StJ}  (24) 


The  production  of  turbulent  energy  becomes 


P  =  ~T  *(*»[ 


dU{i) 

dx(j) 


^U(j)Hj(i)  +  U(l)Hi{l)6i, 


-  3  \pk  +  peffV(i)U(i)}V(l)U(l) 


(25) 


after  the  transformation  to  general  orthogonal  coordinates.  All  of  the  transport 
equations  can  be  written  as 


dxl'  dxl 

V(i)[Pf/(iK’- 

r*o)  *_/i\  ]  +  V(2)[pU(2)ip  -  r 

*(2) 


’ax(l)' 


+  V(3)[/>t'(3)v  -  r  i(3) 


dx(2)1 
dn'  . 


(26) 


&r(3)‘ 


= 


where  all  and  S ^  are  defined  in  Table  1  for  the  appropriate  ip,  and  U,V,W 

correspond  to  U (1),U (2),U (3). 


The  algebraic  finite  difference  equations  can  be  easily  derived  by  applying  the 
Gauss  Divergence  theorem  to  Eq.  26-  The  integration  is  performed  over  the  entire 


13 


control  volume.  Figure  4  depicts  a  general  three-dimensional  control  volume.  The 
variable  tp  is  stored  at  nodes  P,N,S,E,W,T,  and  B.  Applying  the  theorem  yields 


ICC 


pU(l)ip-  rx-^dx{2)dx{Z) 


F y  ~C~dz{l)dx(Z) 


''Apr'' 


pU{2)%p~ 


\\PP 

xi  Uxi  Jx  i 


pU(Z)rP-Tz^dx(l)dx(2)}  +  = 

OZ  ,3 


.3  .,3  . _ 1 


r*+  t*+  t*+ 

/  /  /  S^dx(l)dx(2)dz(3) 

J  X*  J  Z*  J  X* 


where  the  left  side  contains  convection  and  diffusion  and  the  right-hand  side 
contains  the  source  term.  This  form  is  similar  to  that  for  Cartesian  coordinates. 
Note  that  dx(l),  dx(2),  and  dx( 3)  are  physical  lengths  in  the  general  orthogonal 
coordinate  system.  Thus,  the  appropriate  scale  factors  are  required  to  compute  the 
integrals.  Upon  integration  the  resulting  algebraic  finite  difference  equation  is  of 
the  form 

ii* p(An  +  As  +  Ae  +  Aw  +  At  -f  Ab)  =  AnV'n  -f  AsV’s 

(28) 

+Ae<Ae  +  AwV’w  +  AtV’t  +  AbV'b  +  (S^)-pvol 

It  is  well  known  that  the  hybrid  differencing  scheme  produces  false  diffusion 
under  certain  conditions  when  large  control  volume  Peclet  numbers  exist.  False 
diffusion  occurs  if  the  streamlines  of  the  flowfield  are  oblique  to  the  grid  lines  and 
a  nonzero  gradient  of  the  dependent  variable  normal  to  the  flow  exists.  This  false 
diffusion  is  a  truncation  error  in  the  finite  difference  formulation.  The  QUICK 
differencing  scheme  of  Leonard  [19]  generally  reduces  false  diffusion.  Rhode  et  al 
20  have  previously  implemented  QUICK  into  a  two-dimensional  labyrinth  seal 
code.  Substitution  of  the  corresponding  interpolation  functions  into  F,q.  27  ,  or  use 


of  the  hybrid  upwind-central  differencing  scheme,  yields 


VI? 


Ei  A ii'i  +  Su 
E,  Ai  -  Sp 


(29) 


where  i  =E,W,N,S,T,B,EE,WW,NN,SS  and  the  last  four  neighbors  are  not  needed 
if  the  hybrid  scheme  is  used. 

A  system  of  staggered  grids  is  used  to  store  the  variables  of  interest  as  is  done  in 
the  TEACH  code  [21].  The  values  of  pressure,  turbulent  kinetic  energy,  turbulence 
dissipation,  and  enthalpy  are  stored  at  the  intersection  points  of  the  primary  grid, 
whereas  each  velocity  component  is  stored  on  a  separate  grid. 

The  boundary  conditions  play  an  extremely  important  role  in  determining 
the  solution  to  flow  in  the  labyrinth  seal.  Wall  functions  based  on  the  log 
law  of  the  wall  are  used  to  determine  the  appropriate  shear  stress  near  a  wall. 
Axial  velocity,  circumferential  velocity,  turbulence  kinetic  energy,  and  turbulence 
dissipation  profiles  were  prescribed  at  the  inlet.  All  of  these  inlet  profiles  were 
assumed  to  be  uniform  radially.  The  magnitude  of  these  quantities,  except  for 
circumferential  velocity,  was  determined  by  executing  the  program  with  equivalent 
conditions  at  the  inlet  and  exit.  This  case  corresponds  to  a  fully  developed,  i.e. 
streamwise  periodic  cavity.  From  this  computer  run  the  circumferential  variation  of 
inlet  axial  velocity,  turbulence  energy,  and  energy  dissipation,  due  to  eccentricity, 
were  all  determined.  The  inlet  circumferential  velocity  profile  was  also  uniform  in 
the  radial  direction  and  varied  from  case  to  case.  At  the  circumferential  location  of 
largest  clearance,  this  quantity  was  specified  according  to  the  problem  of  interest. 
All  other  circumferential  locations  were  given  an  inlet  value  such  that  the  swirl 
velocity  conserved  mass  in  the  circumferential  direction.  The  radial  velocity  was 
set  to  zero  at  the  inlet  for  all  circumferential  locations.  Pressure  was  prescribed  at 
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one  nodal  location  and  its  value  held  constant  during  the  numerical  procedure.  At 
this  location  pressure  had  a  value  of  3  atmospheres.  All  other  pressure  nodes  were 
calculated  on  a  relative  basis  with  respect  to  this  node.  The  stagnation  enthalpy  at 
the  inlet  was  calculated  using  the  inlet  velocities  and  temperature. 

Convection  is  dominant  over  diffusion  at  the  exit,  and  the  Hybrid  up¬ 
wind/central  differencing  scheme  applied  there  requires  no  downstream  boundary 
values.  The  axial  velocity  at  the  exit  was  set  equal  to  the  axial  velocity  at  the  up¬ 
stream  node  plus  a  uniformly  distributed  percentage  of  the  axial  velocity  increment 
required  to  conserve  mass  globally. 


RESULTS  AND  DISCUSSION 

The  first  cavity  of  a  multi-cavity  seal  is  considered  throughout  this  work. 
Figure  5  depicts  the  cavity  and  shows  the  relevant  dimensions.  These  dimensions 
were  chosen  because  they  constitute  a  generic  cavity  which  is  somewhat  similar 
to  those  used  in  high  performance  turbomachines.  Further,  this  geometry  can  be 
modelled  with  a  uniform  grid  within  the  cavity.  Due  to  the  coordinate  system 
arrangement  the  grid  must  be  nonuniform  in  the  gap  region.  Minimizing  the 
grid  nonuniformity  as  much  as  possible  was  desireable  in  aiding  convergence.  The 
eccentricity  of  the  rotor  was  three  percent  of  the  clearance  c  in  this  analysis. 

Computational  Characteristics 

A  total  of  six  converged  solutions  were  obtained.  In  all  six  cases  the  rotational 
speed  of  the  rotor  was  5,000  cpm.  The  working  fluid  was  air  at  3.0  atm.  and 
294  K  at  the  inlet.  The  axial  Reynolds  number  Re  =  2U%nc/v  in  all  cases  was 
19,200,  which  corresponds  to  a  Mach  number  near  0.2.  The  inlet  conditions  in 
all  six  cases  were  identical  except  for  swirl  velocity.  Three  different  inlet  swirl 
cases  were  investigated.  These  were  approximately  30,60,  and  90  percent  of  the 
rotational  speed  of  the  cavity,  which  correspond  to  tangential  Reynolds  numbers 
Re#  =  2Winc/v  of  1,315,  2,631,  and  3,946.  Each  swirl  case  was  considered  with  no 
whirl  as  well  as  half-speed  forward  whirl,  i.e.  2,500  cpm. 

A  grid  dependence  study  revealed  that  a  grid  of  22x22x17  nodes  in  the  x-,  r-, 
and  ^-directions,  respectively  was  sufficient  to  simulate  the  flowfield  for  this  seal 
problem.  Only  the  number  of  grid  points  in  the  circumferential  direction  was  varied 
in  the  present  grid  independence  study.  In  order  to  conserve  computer  resources, 
the  number  of  grid  nodes  used  in  the  axial  and  radial  directions,  22x22,  was  chosen 
based  on  extensive  previous  grid  dependence  testing  [22]  at  the  same  Mach  number 
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using  a  concentric-rotor  version  of  the  computer  program.  Grids  of  22x22x11, 
22x22x17,  and  22x22x25  were  all  used  to  solve  the  eccentric  cavity  whirling  at  half 
speed.  All  parameters  used  in  the  grid  dependence  study  were  similar  to  those  used 
in  the  final  computations.  The  22x22x17  grid  was  required  to  determine  the  grid 
independent  value  for  pressure  drop  across  the  cavity.  The  grid-independent  results 
are  presented  in  Figure  6. 

Results 

The  distribution  of  the  tangential  force  component  within  the  cavity  can  be 
observed  in  Figure  7.  This  plot  clearly  shows  that  the  downstream  tooth  contributes 
heavily  to  the  total  tangential  force.  Further,  the  effect  of  inlet  swirl  strength  on 
this  force  component  is  only  observed  at  the  downstream  tooth.  This  is  attributed 
to  the  fact  that  the  circumferential  flow  is  developing  within  the  cavity  so  that  there 
is  a  different  pressure  distribution  at  the  downstream  tooth. 

For  the  cases  considered  here,  the  radial  force  generally  decreases  slightly  with 
increasing  swirl.  Note  that  this  force  component  is  not  restoring  the  rotor  to  its 
original  position,  but  is  actually  pushing  it  to  larger  eccentricities.  These  same 
general  trends  were  measured  by  Leong  [23]  at  small  eccentricities.  As  with  the 
tangential  force,  the  flow  over  the  downstream  tooth  significantly  influences  the  net 
radial  force  as  seen  in  Figure  8.  Observe  that  swirl  has  little  effect  on  the  radial 
force  distribution  and  that  this  force  is  smallest  in  the  middle  of  the  cavity,  as  for 
the  tangential  force. 

It  was  determined  that  the  net  shear  stress  contribution  to  the  total  net  radial 
force  is  less  than  one  percent,  as  is  the  case  with  the  net  tangential  force.  As 
expected,  the  shear  contribution  to  the  radial  force  was  a  minimum  when  Mr*n  =  0.6 
because  this  is  near  the  asymptotic  value  in  which  case  there  was  negligible  Tre  stress 


on  the  south  wall. 

As  observed  earlier,  the  pressure  force  virtually  determined  the  net  forces  acting 
on  a  seal.  The  variation  of  pressure  (relative  to  a  reference  value  given  below) 
in  the  circumferential  direction  at  three  locations  along  the  rotor  are  plotted  in 
Figures  9  through  11.  The  locations  are  on  the  upstream  tooth,  in  the  middle 
of  the  cavity,  and  on  the  downstream  tooth.  The  reference  pressures  used  to 
nondimensionalize  the  pressure  in  the  polar  plots,  are  306,000  Pa,  319,400  Pa,  and 
321,800  Pa  respectively. 

As  seen  in  Figure  9,  swirl  has  little  effect  on  pressure  at  the  upstream  tooth. 
The  distribution  of  pressure  is  somewhat  symmetrical  about  the  axis  through  9  =180 
degrees.  This  symmetry  indicates  that  the  tangential  force  component  is  smaller 
than  the  radial  component  at  the  upstream  tooth. 

In  Figure  10  the  pressure  values  are  taken  at  the  middle  of  the  south  wall  of  the 
cavity.  Symmetry  about  the  same  axis  is  also  observed,  and  swirl  had  little  effect 
on  this.  At  this  mid-cavity  location  the  tangential  and  radial  forces  were  earlier 
seen  to  be  at  a  minimum.  Thus  this  circumferential  pressure  distribution  is  the 
closest  to  being  perfectly  symmetrical  about  the  axis  through  9  =180  degrees  and 
that  through  9  =90  degrees. 

In  Figure  11  the  pressure  distribution  on  the  downstream  tooth  is  given. 
Observe  that  the  pressure  distribution  is  not  symmetric  at  all.  This  shows  the  large 
pressure  forces  on  the  downstream  tooth  seen  earlier  in  Figures  7  and  8.  Higher 
inlet  swirl  values  increased  this  lack  of  symmetry,  which  caused  an  increase  in  the 
tangential  force. 

Although  turbomachines  frequently  whirl  at  half-speed,  the  case  of  no  whirl 
has  also  been  studied  previously  [11.12,13].  The  F*  distribution  along  the  cavity 


for  the  no-whirl  case  was  plotted  in  Figure  12.  These  results  are  very  similar  to 
those  for  half-speed  whirl  shown  in  Figure  7.  This  was  expected  since  for  such  small 
eccentricities  the  boundary  conditions  are  nearly  the  same.  Figure  13  is  the  no¬ 
whirl  counterpart  of  Figure  8  for  F* .  Again,  quantitatively  as  well  as  qualitatively 
the  results  are  very  similar  to  the  whirling  case. 

The  Shear  Stress 

The  shear  stress  distributions  along  the  cavity  walls  provide  previously  unavail¬ 
able  insight  into  the  basic  flowfield.  Although  they  were  not  found  to  contribute 
significantly  to  the  force  components,  estimates  of  these  are  useful  in  the  develop¬ 
ment  of  simpler  models.  Shear  stresses  are  plotted  along  the  wall  and  in  the  free 
shear  layer  for  a  specific  theta  value.  The  circumferential  variation  for  a  specific 
wall  location  is  also  given.  Thus  Figures  14  through  18  can  be  used  to  construct 
the  entire  three-dimensional  shear  stress  distribution  along  the  cavity  walls  and  the 
free  shear  layer.  Figure  14  shows  the  points  labeled  E,F,G,  and  H  which  are  used 
to  denote  a  specific  location  in  the  cavity.  For  example,  wall  FG  is  the  south  wall 
of  the  cavity,  and  EH  is  the  free  shear  layer. 

The  variation  of  shear  stress  along  the  rotor  surface  at  0=0  degrees  for  all 
three  inlet  swirl  cases  at  half-speed  whirl  can  be  seen  in  Figure  15.  Although  the 
rrz  component  does  not  vary  with  inlet  swirl,  ttq  varies  substantially  and  is  largest 
for  the  low  inlet  swirl  case.  This  is  because  the  dw/dr  term  in  the  expression  for 
tT0  is  largest  for  low  inlet  swirl,  as  the  shaft  tangential  velocity  is  much  higher  than 
that  of  the  fluid.  For  the  case  of  high  inlet  swirl,  ttq  is  almost  uniform.  Note  that 
the  stresses  are  smallest  in  the  corners  of  the  cavity  where  the  turbulence  energy 
and  thus  the  turbulent  viscosity  is  small. 

In  Figure  16  the  free  shear  layer  stresses  are  shown.  Again,  Trz  does  not  vary 
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with  inlet  swirl,  and  stress  Tr$  is  fairly  uniform  for  the  high  swirl  case.  The  Tr$ 
stress  is  highest  near  the  upstream  tooth,  as  this  is  again  a  region  of  very  large 
dW/dr. 

The  shear  stresses  along  the  west  wall  of  the  cavity  at  0=0  degrees  can  be  seen 
in  Figure  17.  Stress  Txr  does  not  vary  with  inlet  swirl,  and  rxg  is  largest  for  the 
low  inlet  swirl  case  as  expected.  Stress  txq  is  rather  uniform  for  high  swirl,  and 
the  maximum  value  oi  txq  occurs  near  the  shear  layer.  The  shear  stresses  along 
the  east  wall  of  the  cavity  at  0=0  degrees  are  presented  in  Figure  18.  The  basic 
characteristics  are  generally  similar  to  those  in  Figure  17.  That  is  Trx  does  not  vary 
with  inlet  swirl,  and  rxe  is  largest  for  low  inlet  swirl.  The  maximum  of  stress  txo  is 
near  the  shear  layer  due  to  the  high  turbulence  level  there. 
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CONCLUSIONS 

A  computer  program  has  been  developed  for  the  prediction  of  rotordynamic 
forces  on  an  eccentric  whirling  cavity.  Eccentricity,  shaft  speed,  whirl  rate,  and 
various  inlet  conditions  can  be  specified  by  the  user.  This  work  demonstrated  that 
numerical  methods  can  be  used  to  determine  rotordynamic  forces  on  iabryinth  seals, 
as  well  as  the  bulk  flow  analytic  models  used  in  the  past.  The  solution  obtained  from 
the  computer  code  appeared  to  be  reasonable  and  contained  the  expected  flowfield 
features.  Comparison  of  forces  with  Leong  [22]  showed  qualitative  agreement. 
Unfortunately,  force  measurements  for  a  single  teeth-on-rotor  labyrinth  cavity  are 
not  available  for  quantitative  comparison.  However,  multi-cavity  predictions  for 
comparison  with  measurements  are  forthcoming. 

For  the  present  case  of  one  generic  cavity  at  three  percent  eccentricity,  rotating 
at  5,000  r.p.m.,  whirling  at  2,500  r.p.m.,  and  with  an  inlet  axial  velocity  near  Mach 
0.2,  the  following  can  be  asserted: 

1.  An  increase  in  inlet  swirl  velocity  increased  the  tangential  force.  The 
direction  of  the  force  would  accelerate  the  whirl. 

2.  The  largest  contribution  by  far  to  the  total  tangential  force  was  from  the 
pressure  imbalance  on  the  periphery  of  the  downstream  tooth.  It  is  expected 
that  the  tangential  force  can  be  affected  by  altering  the  geometry  of  the  cavity 
at  the  downstream  tooth.  Specifically,  the  thickness  and/or  geometry  of  the  tooth 
periphery  is  apparently  a  very  important  consideration  when  one  is  designing  a  bank 
of  labyrinth  seal  cavities  in  high  performance  turbomachinery  where  rotordynamic 
instabilities  are  a  problem. 

3.  The  radial  force  was  not  restoring  the  shaft,  but  was  pushing  the  shaft  to 
larger  eccentricities. 
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Table  I  Diffusion  Coefficients  and  Source  Terms  in  the  Transport  Equations 
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Fig.  1  Resulting  streamlines  due  to  leakage  flow  in  a  labyrinth  seal 
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Fig.  8  Variation  of  total  radial  force  along  the  rotating 
surface  for  three  inlet  swirl  cases  of  half-speed  whirl 
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Fie.  16  Variation  of  shear  stresses  along  the 
free  shear  layer  at  rt  -  0  decrees  for  cases  A.B,  and  C 
( W •  =0.3. 0.6.  and  0.9  respectively  )  at  half-speed  whirl 
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Fig.  17  Variation  of  shear  stresses  along  the 
upstream  tooth  at  9  =  0  degrees  for  cases  A.B,  and  C 
(Hr'n  =0.3. 0.6.  and  0.9  respectively)  at  half-speed  whirl 
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Fig.  18  Variation  of  .shear  stresses  along  the 
downstream  tooth  at  9  =  0  degrees  for  cases  A,B.  and  C 
(ir.*n  =0.3.0. 6,  and  0.9  respectively)  at  half-speed  whirl 
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ABSTRACT 

A  finite  difference  method  for  determining  rotordynamic  forces  on  an  eccentric 
whirling  labyrinth  cavity  has  been  developed.  A  coordinate  transformation  was 
applied  to  the  turbulent  flow  Navier-Stokes  equations  in  order  to  use  the  modified 
bipolar  coordinate  system.  The  SIMPLER  algorithm  with  QUICK  differencing  and 
the  high  Reynolds  number  k-£  turbulence  model  are  used  to  compute  the  complex 
turbulent  flowfield.  A  circular  whirl  orbit  about  the  geometric  center  of  the  housing 
was  specified  for  simplicity.  For  the  cases  considered,  the  radial  and  tangential 
force  components  increased  and  decreased,  respectively,  with  increasing  inlet  swirl. 
Also,  circumferential  pressure  variations  are  included  for  enhanced  insight  into  the 
flowfield.  Further,  the  circumferential  variation  of  both  shear  stress  components 
along  each  surface  of  the  cavity  are  presented  to  allow  the  developers  of  various 
bulk  flow  models  to  refine  their  stress  modelling. 
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INTRODUCTION 

Self-ex  cited  vibration  remains  a  serious  problem  encountered  in  high  perfor¬ 
mance  turbomachines.  A  complete  quantitative  analysis  of  the  excitation  forces 
causing  such  vibration  is  still  not  complete.  The  excitation  forces  result  from  the 
displacement  of  the  rotor  from  its  geometric  center  in  the  housing.  This  displace¬ 
ment  often  results  from  a  small  mass  imbalance  on  the  rotor.  The  result  of  the 
self-excited  forces  is  usually  subsynchronous  whirl.  Subsynchronous  whirl  (whirling 
at  a  frequency  below  the  rotational  frequency  of  the  shaft)  may  occur  in  the  same 
direction  as  the  shaft  rotation  or  in  the  opposite  direction.  One  source  of  these 
self-excited  forces  is  labyrinth  seals. 

Labyrinth  seals  are  used  to  limit  internal  leakage  from  high  to  adjacent  low 
pressure  regions  in  a  turbomachine.  They  consist  of  a  highly  frictional  flow  path 
between  the  rotor  and  stator.  The  flow  path  converts  pressure  head  into  mean  flow 
kinetic  energy,  which  is  dissipated  by  viscous  losses.  Figure  1  shows  the  streamlines 
within  a  short  bank  of  labyrinth  cavities. 

A  number  of  analytical  investigations  on  the  destablizing  forces  generated  by 
the  leakage  flow  through  labyrinth  seals  have  been  performed.  Many  of  the  analyti¬ 
cal  methods  employ  bulk  flow  models  and  use  numerous  empirical  correlations.  The 
problem  was  first  modelled  by  Alford  [1],  who  neglected  circumferential  flow,  result¬ 
ing  in  no  predicted  force.  The  calculations  of  Kostyuk  [2]  neglected  the  change  of 
chamber  area  due  to  rotor  eccentricity  and  contradicted  various  force  measurements. 
Iwatsubo  [3,4]  extended  Kostyuk’s  model  to  include  this  chamber  area  variation. 
Childs  and  Scharrer  [5]  included  this  variation  and  developed  improved  shear  stress 
relationships.  They  obtained  generally  improved  predictions.  These  investigators 
also  presented  [6]  the  first  complete  set  of  separate  experimental  values  for  the 


stiffness  and  damping  coefficients. 

Only  very  recently  has  Computational  Fluid  Dynamics  CFD  been  utilized 
to  predict  the  forces  on  labyrinth  seals.  Dietzen  and  Nordman  [7]  used  it  in  a 
hybrid  perturbation/CFD  model  to  solve  the  Navier-Stokes  equations  for  turbulent 
flow.  Wyssman  [8]  was  apparently  the  first  to  present  pressure  distributions  for 
an  eccentric-rotor  labyrinth  cavity  from  a  three-dimensional  computation  of  the 
turbulent  Navier-Stokes  equations.  These  along  with  velocity  distributions  given 
by  Wyssman,  et  al  [9]  were  used  to  develop  a  special  two-control-volume  bulk-flow 
model.  Tam,  et  al  [10]  used  a  three-dimensional  CFD  code  to  compute  the  forces 
in  non-labyrinth,  i.e.  annular  seals. 

NUMERICAL  MODEL 

In  this  study  the  whirl  orbit  has  been  idealized  as  circular  motion,  the  center 
of  which  coincides  with  that  of  the  housing.  A  single  cavity  of  a  multi-cavity  seal 
of  teeth-on-rotor  design  is  analyzed.  Since  the  tooth  perimeter  is  concentric  with 
the  base  of  the  seal  cavity,  cylindrical  coordinates  are  the  natural  choice  within 
the  seal  cavity.  However,  in  the  clearance  region  between  the  teeth  and  the  stator, 
bipolar  coordinates  are  used  because  of  the  eccentricity  of  the  shaft.  Specifically, 
the  modified  bipolar  coordinate  system  was  selected  as  it  becomes  the  cylindrical 
coordinate  system  in  the  limit  as  eccentricity  approaches  zero.  Figure  2  shows 
such  a  finite  difference  grid  within  the  clearance  region.  Reference  [11]  gives  details 
concerning  the  modified  bipolar  coordinate  system. 

Neglecting  body  forces,  the  steady-state  continuity  and  momentum  equations 


where  p,  U,,  and  p  are  the  time-averaged  density,  velocity  and  pressure.  The 
Reynolds  stresses  are  determined  by  the  high  Reynolds  number  k  —  e  turbulence 
model. 


The  equations  were  transformed  into  relations  involving  general  orthogonal 
coordinates,  using  Pope’s  [12]  method.  With  his  transformation  procedure  the 
equations  retain  their  original  form  and  simplicity  as  much  as  possible. 

Consider  a  general  orthogonal  coordinate  system  where  the  orthogonal  coordi¬ 
nates  are  denoted  by  x‘.  Distances  in  this  system  can  be  related  to  the  Cartesian 
system  by  the  scale  factors  hi, 

(ds)2  =  (h.dx*)2  =  [dx(i)]2  (3) 


where  the  scale  factors  h,  are  excluded  from  the  summation  convention.  Thus, 
the  physical  displacements  along  a  coordinate  line  x*  in  the  general  orthogonal 
coordinate  system  are  dx(i).  The  scale  factors  are  determined  by 


2  _  dxi  dxi 
1  ~  dxl  dxl 


(4) 


where  x\  are  the  coordinates  in  the  Cartesian  coordinate  system.  Pope’s  diver¬ 
gence  operator  V(i)  and  coordinate  variation  term  Ht(j )  needed  to  transform  the 
transport  equations  into  general  orthogonal  coordinates  are 
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Here  /i,  represent  scale  factors  and  jh|  is  the  product  of  the  scale  factors.  After 
performing  the  transformation,  the  mass  and  momentum  equations  become 

V(i)[pU(i)}  =  0  (7) 


V(i)[pf/(i)t7(i)  +  r*(tj)]  =  --f-  +/f,(i)[pf/(0l/(0  +  r*(«)] 

dx(j)  (8) 

-Bi(i)\pU(ipU)  +  rm{ij)) 

The  isotropic  component  of  stress  has  been  added  to  the  pressure  giving 

V  =  P  +  ^pk+  |/i«//V(*)l7(t)  (9) 

where  r*  has  the  anisotropic  stress  as 

T-(.j)  =  -  V(i)n,(i)  +  2um,mi)  (10) 

All  of  the  transport  equations  can  be  written  as 

V(1)1PU(1W  -  r,(1)^L]  +  V(2)[pt/(2)V>  -  rl(2)-^-] 

o\  (11) 

+  V(3)[pU(3)V>-r,(3)^]  =  S* 

where  all  and  S $  are  defined  in  reference  [11]  for  the  appropriate  ip,  and 

U,  V ,  W  correspond  to  U(1),U(2),U(Z). 

The  algebraic  finite  difference  equations  can  be  easily  derived  by  applying  the 
Gauss  Divergence  theorem  to  Eq.  11.  Further  details  are  available  in  reference  [11]. 
Figure  2  shows  the  modified  bipolar  coordinate  portion  of  the  grid.  Observe  that 
di(l),dx(2),  and  dx( 3)  are  physical  lengths  in  the  general  orthogonal  coordinate 
system.  Upon  integration  the  resulting  algebraic  finite  difference  equation  is  of  the 


i/’p(Afj  +  As  +  Af?  +  Aw  +  Aj  +  Ag)  =  AnV’n  +  AsV’s 
+Ae«/,e  +  AwV’w  +  At*/’t  +  AbV’b  +  (S^)pvol 
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The  hybrid  upwind/central  differencing  scheme  yields  false  diffusion  under 
certain  conditions  when  large  control  volume  Peclet  numbers  occur.  False  diffusion 
error  arises  if  the  streamlines  are  oblique  to  the  grid  lines  and  a  nonzero  gradient  of 
the  dependent  variable  normal  to  the  flow  exists.  The  QUICK  differencing  scheme 
of  Leonard  [13]  generally  reduces  false  diffusion.  Rhode  et  al  [14]  have  previously 
implemented  QUICK  into  an  axisymmetric  computer  code  from  which  this  program 
was  developed.  The  three-point  interpolation  expression  for  the  west  face  of  the 
control  volume  on  a  uniform  grid  is  shown  in  Figure  3  and  is  calculated  by 


V’w  = -(Vr  +  V*w)  ~  g(V,E  -  2 V'P+V’w)  (13) 

if  Uw  is  negative  and 

V’w  =  —  ( +  tAw)  -  ~(t^p  -  2t/>w  +  V’ww)  (14) 

1  O 


if  Uw  is  positive.  Tv  interpolation  functions  were  modified  for  any  non-uniform 
grid.  Substitution  of  these  interpolation  functions  into  Eq.  12  ,  or  use  of  the  hybrid 
differencing  scheme  yields 


V>p  = 


S.  A.V’i  +  Su 
Si  Ai  -  Sp 


(15) 


where  i  =E,W,N,S,T,B,EE,WW,NN,SS  and  the  last  four  neighbors  are  not  needed 
in  th  *  hybrid  scheme. 

As  is  done  in  the  TEACH  code  [15],  the  variables  were  stored  on  a  system  of 
four  staggered  grids.  The  values  of  pressure,  turbulent  kinetic  energy,  turbulence 
dissipation,  and  enthalpy  were  stored  at  the  intersection  points  of  the  primary  grid, 
whereas  each  velocity  component  was  stored  on  a  separate  grid. 


Wall  functions  were  used  to  determine  the  appropriate  shear  stress  near 
a  wall.  A  radially  unifrom  profile  of  axial  velocity,  circumferential  velocity. 


turbulence  kinetic  energy,  and  turbulence  dissipation  was  prescribed  at  the  inlet. 
The  circumferential  variation  of  each  was  determined  from  a  preliminary  computer 
rim.  For  that  run,  the  inlet  value  of  each  variable,  except  pressure  and  swirl  velocity, 
was  set  equal  to  the  corresponding  outlet  value  from  the  previous  iteration.  This 
case  corresponds  to  a  fully  developed,  i.e.  st  re  am  wise  periodic  cavity.  The  inlet 
circumferential  velocity  profile  was  also  uniform  in  the  radial  direction.  At  the 
circumferential  location  of  largest  clearance,  this  quantity  was  specified  according 
to  the  problem  of  interest.  All  other  circumferential  locations  were  given  an  inlet 
value  such  that  the  swirl  velocity  conserved  mass  in  the  circumferential  direction. 
Pressure  was  prescribed  at  one  point  in  the  domain. 

The  solution  strategy  employed  is  the  well  documented  SIMPLER  procedure. 
A  line  solver  is  used  to  solve  for  each  variable  of  interest.  This  line  solver  is  applied 
in  the  ^-direction  in  the  cylindrical  coordinate  system  and  in  the  ^-direction  in  the 
modified  bipolar  coordinate  system.  The  re-entrant  or  cyclic  boundary  condition 
is  utilized.  The  assembled  equations  are  written  in  the  form 

-AtV’t  +  (Ap  -  Sp)^p  -  AbV’b  =  Su  +  ^  A,t/>,  (16) 

i 

where  i  =N,S,E,W,NN,SS,EE,WW.  Application  of  Eq.  16  to  each  grid  circle  leads 
to  a  matrix  that  is  solved  using  a  CTDMA  algorithm.  Each  grid  circle  throughout 
the  domain  is  solved  in  this  manner. 


DISCUSSION 

Figure  4  shows  the  cavity  and  the  relevant  dimensions.  Minimizing  the  grid 
nonunifonnity  was  desireable  in  aiding  convergence.  The  eccentricity  of  the  rotor 
was  three  percent  of  the  clearance  c  in  this  analysis. 

TheFlowfieldsConsidered 

Three  different  inlet  swirl  cases  were  investigated:  30,60,  and  90  percent  of 
the  rotational  speed  of  the  cavity.  The  rotational  speed  of  the  rotor  was  5,000 
cpm  and  the  working  fluid  was  air  at  3.0  atm.  and  294  K  at  the  inlet.  The  axial 
Reynolds  number  Re*  =  2Uinc/u  in  all  cases  was  19,200,  which  corresponds  to  a 
Mach  number  near  0.2.  The  inlet  conditions  were  identical  for  each  case  except  for 
swirl  velocity.  Each  swirl  case  was  considered  with  no  whirl  as  well  as  half-speed 
forward  whirl  at  2,500  cpm. 

A  grid  of  22x22x17  nodes  in  the  x-,  r-,  and  ^-directions,  was  found  in  a 
grid-dependence  study  to  be  sufficient  to  simulate  the  flowfield  for  this  problem. 
A  considerably  finer  grid  is  required,  however,  for  higher  Mach  number  cases. 
In  order  to  conserve  computer  resources,  the  number  of  grid  nodes  used  in  the 
axial  and  radial  directions,  22x22,  was  determined  in  conjuction  with  previous  grid 
dependence  testing  [16]  at  the  same  Mach  number. 

For  the  present  computer  runs,  it  was  found  that  the  radial  component  of 
the  net  force  was  much  larger  than  the  tangential  component.  Inasmuch  as 
the  tangential  component  is  a  very  small  value  here,  it  meandered  somewhat 
(approximately  *  8  percent)  with  further  iterations  beyond  the  point  of  convergence 
for  other  quantities. 

Results 

The  rotor  in  actual  turbomachines  has  frequently  become  unstable  with  the 


shaft  whirling  at  a  frequency  approximately  half  of  the  rotational  frequency.  Results 
for  this  important  case  of  half-speed  whirl  are  considered  here.  The  net  tangential 
force  acting  on  the  eccentric  whirling  shaft  is  of  considerable  interest  because  it  is 
the  driving  force  for  the  whirl  phenomenon.  In  Figure  5  the  effect  of  inlet  swirl 
strength  on  the  net  tangential  force  can  be  seen.  The  total  net  tangential  force 
contains  contributions  from  both  pressure  and  shear  forces.  As  expected  from  the 
literature,  increasing  swirl  increases  the  driving  force  behind  the  rotordynamic  whirl. 

The  contribution  of  shear  forces  to  the  total  net  tangential  force  is  seen  in 
Figure  6.  From  the  quantitative  result  in  Figures  5  and  6  it  is  obvious  that  the 
shear  component  of  the  net  tangential  force  is  less  than  one  percent.  The  other  99 
percent  is  from  the  net  pressure  force  due  to  the  circumferential  pressure  distribution 
around  the  cavity.  Observe  that  the  south  wall(cavity  base)  contribution  is  less  than 
zero  and  the  others  are  greater  than  zero.  This  gives  a  cancellation  and  a  small 
total  shear  force.  Also,  note  that  the  south  wall  shear  was  the  only  one  significantly 
influenced  by  a  change  in  inlet  swirl. 

While  a  positive  tangential  force  is  a  whirl  driving  force,  a  negative  radial  force 
is  a  restoring  force.  The  effect  of  inlet  swirl  on  the  total  net  radial  force  can  be  seen 
in  Figure  7.  The  radial  force  generally  decreases  slightly  with  increasing  swirl.  It  is 
about  four  times  larger  than  the  tangential  force.  Note  that  the  force  is  not  restoring 
the  rotor  to  its  original  position,  but  is  actually  pushing  it  to  larger  eccentricities. 
These  same  general  trends  were  measured  by  Leong  [17]  at  small  eccentricities. 

It  was  determined  from  Figures  7  and  8  that  the  net  shear  stress  contribution 
to  the  total  net  radial  force  is  less  than  one  percent,  as  is  the  case  with  the  net 
tangential  force.  Figure  8  also  shows  that  all  radial  shear  force  components  were 
invariant  with  swirl,  except  at  the  south  wall.  As  expected,  the  shear  contribution 


to  the  radial  force  there  was  a  minimum  when  W*n  =  0.6  because  this  is  near  the 
asymptotic  value  in  which  case  there  was  negligible  Tr$  stress  on  the  south  wall. 


Since  the  pressure  force  virtually  determines  the  net  forces  acting  on  a  seal,  the 
variation  of  pressure  (relative  to  a  reference  of  306,000  Pa)  in  the  circumferential 
direction  at  two  locations  along  the  rotor  are  plotted.  In  Figure  9  the  pressure 
values  are  taken  at  the  middle  of  the  south  wall  of  the  cavity.  Symmetry  about 
the  axis  through  0=180  degrees  is  observed  indicating  a  very  small  contribution 
to  Ft,  and  swirl  had  little  effect  on  this.  In  Figure  10  the  pressure  distribution 
on  the  downstream  tooth  is  given.  Observe  that  here  the  pressure  distribution  is 
not  symmetric  at  all  through  this  axis,  contributing  substantially  to  Ft-  Higher 
inlet  swirl  values  produced  this  lack  of  symmetry,  which  caused  an  increase  in  the 
tangential  force. 

The  case  of  no  whirl  has  also  been  studied  previously.  The  current  results  for 
the  non-whirling  case  are  very  similar  to  those  shown  for  half-speed  whirl.  This 
was  expected  since  for  such  small  eccentricities  the  boundary  conditions  are  nearly 
the  same. 

The  development  of  circumferential  velocity  within  a  bank  of  labyrinth  cavities 
has  been  of  interest  for  some  time  because,  as  shown  above,  swirl  has  a  substantial 
influence  on  the  forces  generated  by  an  eccentric  seal.  The  circumferential  variation 
of  swirl  at  the  inlet  and  exit  for  the  case  of  lowest  inlet  swirl  is  given  in  Figure  11. 
The  swirl  increased  in  this  case  from  the  inlet  to  exit  by  about  seven  percent.  This 
was  expected  since  the  inlet  swirl  was  below  the  fully  developed  (i.e.  asymptotic) 
value.  Notice  that  the  peak  of  the  distribution  shifted  circumferentially  nearly  25 
degrees.  This  is  attributed  to  the  eccentricity  and  swirl  development.  For  the 
intermediate  inlet  swirl  case  of  W*  =  0.6,  the  swirl  velocity  increased  by  four 


15 


percent.  Since  this  value  of  inlet  swirl  is  closer  to  the  asymptotic  limit,  there  was 


less  of  an  increase.  The  peak  value  shifted  only  about  12  degrees,  in  this  instance. 


For  the  largest  inlet  swirl  case,  the  swirl  changed  little  in  magnitude  and  shifted 


about  12  degrees. 


The  development  of  axial  velocity  through  the  cavity  was  also  examined.  This 


development  for  all  three  inlet  swirl  cases  is  presented  in  Figure  12.  In  each  case 


the  axial  velocity  inlet  profile  was  exactly  identical,  and  all  three  cases  produced 


the  same  exit  axial  velocity.  The  exit  axial  velocity  was  almost  symmetric  about 


an  axis  through  9  =180  degrees.  For  the  present  inlet  Mach  number,  the  velocity 


increased  by  about  two  percent  from  inlet  to  exit. 


The  velocity,  pressure,  turbulence  energy  and  turbulence  dissipation  distribu¬ 


tions  give  additional  insight  into  the  flowfield  within  the  cavity  which  is  useful  to 


the  developers  of  of  simpler  models.  Reference  [18]  contains  profile  plots  of  these 


quantities.  These  plots  are  for  the  intermediate  inlet  swirl  case  with  half-speed 


whirl.  The  velocity  distributions  within  the  cavity  were  quite  similar  to  those  pro¬ 


duced  by  Rhode  and  Sobolik  [14]  for  a  concentric-rotor  seal.  The  velocity  did  not 


vary  greatly  circumferentially.  This  can  be  expected  since  the  eccentricity  was  only 


three  percent.  The  axial  velocity  component  shows  significantly  more  circumferen¬ 


tial  variation  than  the  other  two  velocity  components. 


Although  they  do  not  contribute  significantly  to  the  force  components,  shear 


stress  distributions  are  needed  in  the  development  of  simpler  models.  Shear  stresses 


are  plotted  along  each  wall  as  well  as  the  free  shear  layer  for  0=0  in  Reference 


[111.  The  corresponding  circumferential  variation  of  these  quantities  at  a  location 


midpoint  along  each  surface  is  given  here.  Thus  the  two  sets  of  figures  together  can 


be  used  to  approximate  the  entire  three-dimensional  shear  stress  distribution  along 


the  cavity  walls  and  the  free  shear  layer. 

The  circumferential  variation  of  the  above  shear  stresses  for  a  location  midway 
along  each  respective  plane  is  shown  in  Figures  13  through  16.  The  stresses  in 
Figure  13  are  for  the  midpoint  along  the  south  wall  of  the  cavity.  Stress  rri  is 
independent  of  inlet  swirl  and  exhibits  a  minimum  at  0=180  degrees  due  to  the 
lower  axial  velocity  at  that  0-location  (most  of  the  leakage  flow  is  near  0=0).  Stress 
rr0  was  fairly  invariant  with  0,  but  swirl  affected  the  magnitude  of  Tr#  as  found  in 
the  previous  figures. 

The  shear  stresses  midway  along  the  free  shear  layer  are  plotted  in  Figure 
14.  As  with  the  corresponding  south  wall  stresses  in  Figure  13,  component  rrx  is 
independent  of  the  inlet  swirl  and  has  a  minimum  at  0=180  degrees.  Stress  Trg  was 
very  sinusoidal  and  decreasing  inlet  swirl  increased  the  magnitude  of  the  stress,  as 
expected. 

The  circumferential  variation  of  east  and  west  wall  midpoint  shear  stresses  can 
be  observed  in  Figures  15  and  16.  For  both  walls  tti  is  again  independent  of  inlet 
swirl,  and  exhibits  a  minimum  at  0=180  degrees.  Stress  txq  is  nearly  uniform  with 
theta  and  decreases  with  increasing  inlet  swirl. 

CONCLUSIONS 

It  has  been  shown  that  rotordynamic  forces  on  an  eccentric  whirling  labyrinth 
seal  can  be  calculated  using  this  numerical  approach.  A  qualitative  comparison  of 
forces  with  Leong  [17]  gave  good  agreement.  Although  measurements  of  the  rotor- 
dynamic  forces  on  a  single  theeth-on-rotor  cavity  are  not  available  for  quantitative 
comparisons,  numerical  simulations  of  multi-cavity  domains  in  the  near  future  will 
allow  this. 

The  current  computational  results  serve  to  provide  the  seal  designer  wdth 


more  complete  information  about  pressure  circumferential  variation  within  a  cav- 
ity.  Specifically,  the  included  results  quantify  the  dependence  of  the  tangential 
and  radial  force  components  on  the  level  of  inlet  swirl.  For  the  cases  considered 
here,  an  increase  in  inlet  swirl  re-distributed  the  pressure  field  to  give  an  increased 
tangential  force  and  a  decreased  radial  force.  Further,  the  circumferential  pres¬ 
sure  distributions  presented  reveal  details  concerning  the  large  variation  of  those 
force  components  with  axial  location  from  the  mid-cavity  location  to  that  of  the 
downstream  tooth. 

It  was  further  found  that,  while  there  is  definitely  a  net  fluid  shear  force  acting 
on  the  rotor,  the  combined  effect  of  all  stresses  tend  to  cancel,  giving  less  than  a  one 
percent  contribution  to  each  force  component.  The  circumferential  variation  of  both 
wall  and  "free”  shear  stresses  is  presented  in  order  to  allow  the  developers  of  simpler 
flow  models  to  refine  their  shear  stress  modeling.  The  same  stress  components  are 
quantified  in  reference  [11]  as  a  function  of  position  along  each  surface  at  a  fixed 
0-location  so  that  an  approximate  three-dimensional  stress  field  can  be  constructed 
from  the  two  sets  of  stress  data  together. 
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Fig.  6  Effect  of  inlet  swirl  on  the  shear  stress  contribution  to 
the  total  tangential  force  on  the  cavity  walls  for  half-speed  whirl 
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Fig-  8  Effect  of  inlet  swirl  on  the  shear  stress  contribution 
to  the  total  radial  force  on  the  cavity  walls  for  half-speed  whirl 
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Fig-  9  Variation  of  relative  pressure  around  the  periphery  of  the 
rotor  at  mid-cavity  for  three  inlet  swirl  cases  at  half-speed  whirl 
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Fig.  11  Circumferential  variation  of  spatially 
averaged  swirl  velocity  component  at  the  cavity 
inlet  and  outlet  for  W*n  =  0.3  and  half-speed  whirl 
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Fig.  12  Circumferential  variation  of  spatially  averaged 
axial  velocity  component  at  the  cavity  inlet  and  outlet  for  all 
three  inlet  swirl  cases  (Wt*n  =  0.3, 0.6,  and  0.9)  at  half-speed  whirl 
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Fig.  13  Circumferential  variation  of  rotor  surface 
shear  stresses  at  mid-cavity  for  cases  A,B,  and  C 
(HA*n  =0.3, 0.6,  and  0.9  respectively)  at  half-speed  whirl 
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Fig.  14  Circumferential  variation  of  free  shear  layer 
shear  stresses  at  mid-cavity  for  cases  A,B)  and  C 
{W*n  =0.3,0. 6,  and  0.9  respectively)  at  half-speed  whirl 
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Fig.  15  Circumferential  variation  of  tooth  shear  stresses 
at  the  midpoint  along  the  downstream  tooth  for  cases  A,B  and  C 
{W*n  =0.3, 0.6,  and  0.9  respectively)  at  half-speed  whirl 
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Fig.  16  Circumferential  variation  of  shear  stresses 
at  the  midpoint  along  the  upstream  tooth  for  cases  A,B,  and  C 
(lFt*n  =  0.3, 0.6,  and  0.9  respectively)  at  half-speed  whirl 
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point  coefficient 
neighbor  coefficient 
leakage  flow  area 

area  of  east  face  of  control  volume 
area  of  north  face  of  control  volume 
area  of  west  face  of  control  volume 
clearance 

diffusion  coefficient 

turbulence  kinetic  energy  dissipation  rate 

diffusion  coefficient 

ratio  of  specific  heats 

turbulence  kinetic  energy 

turbulence  length  scale 

leakage  mass  flow  rate 

Mach  number 

absolute  viscosity 

effective  viscosity 

turbulent  viscosity 

static  pressure 

pressure  correction 

Peclet  number 

bulk  averaged  pressure 

radial  coordinate 

outer  stator  radius 

Reynolds  number 

density 

bulk  averaged  density 

source  term 

general  source  term 

Taylor  number 

axial  velocity  component 

bulk  averaged  axial  velocity  component 

radial  velocity  component 

specific  volume 

swirl  velocity  component 

axial  coordinate 
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Subscripts 


cavity  exit,  east  face  of  control  volume 
eastern  neighbor  control  volume 
eastern  neighbor  of  eastern  neighbor 
cavity  inlet 

north  face  of  control  volume 

northern  neighbor  control  volume 

northern  neighbor  of  northern  neighbor 

center  control  volume 

center  of  control  volume 

general  flow  variable 

radial 

south  face  of  control  volume 
southern  neighbor  control  volume 
southern  neighbor  of  southern  neighbor 
shaft 

circumferential 

west  face  of  control  volume 

western  neighbor  control  volume 

western  neighbor  of  the  western  neighbor 

axial 


Superscripts 


nondimensionalized,  uncorrected  value 
correction  value 
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ABSTRACT 

A  recently  developed  version  of  a  swirl-flow  finite  difference  computer  program 
was  improved  and  employed  in  predicting  the  compressible  flow  of  air  through 
labyrinth  seeds.  The  substantial  effect  of  inlet  leakage  Mach  number  on  grid 
sensitivity  of  the  solution  is  investigated.  Further,  cavity-by-cavity  development  of 
the  flowfield  is  computed  and  the  distribution  of  various  field  variables  are  presented. 
Results  are  for  straight-through  seals  of  both  teeth-on-stator  and  teeth-on-rotor 
types.  The  teeth-on-rotor  seal  gives  less  leakage  than  the  equivalent  teeth-on-stator 
design.  However,  it  exhibits  a  greater  tendency  to  generate  self  excited  rotordynamic 
forces  due  to  higher  swirl  velocities.  Also,  previously  unavailable  predictions  of 
swirl  velocity  development  are  provided  for  the  refinement  of  simple  models  for 
rotordynamic  forces. 
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INTRODUCTION 

Labyrinth  seals  axe  most  commonly  used  in  rotating  machinery  such  as  pumps, 
compressors  and  turbines.  C.  A.  Parsons  [1]  apparently  was  the  first  to  introduce 
the  labyrinth  seal  in  his  development  of  the  steam  turbine  near  the  turn  of  the 
century.  Figure  1  shows  a  straight-through,  teeth-on-rotor  labyrinth  seal. 

The  sealing  objective  is  to  present  a  highly  frictional  flow  path  between  high- 
and  low-pressure  regions  by  means  of  a  series  of  non-contacting  restrictors  and 
separating  chambers.  When  the  fluid  enters  a  seal,  it  flows  through  a  small 
constriction  at  the  first  tooth  and  part  of  the  pressure  head  is  converted  into  mean 
flow  kinetic  energy.  Seals  are  designed  so  that  a  large  portion  of  this  kinetic  energy 
is  lost  via  turbulence  dissipation  in  the  chamber  immediately  downstream.  Other 
constrictions  and  chambers  follow  downstream,  w’here  the  process  is  repeated. 

Numerous  investigators  have  proposed  empirically  based  relations  utilizing 
characteristics  of  the  overall  flowfield  for  estimating  the  leakage  rate.  Experimental 
data  such  as  total  pressure  drop  has  been  recorded  and  used  to  develop  these 
relationships.  Leakage  has  been  expressed  as  a  function  of  overall  pressure  drop, 
friction  factor,  seal  clearance,  tooth  thickness,  cavity  width,  shaft  speed  and  number 
of  teeth.  Resulting  predictions  from  these  formulas  are  successful  when  applied  to 
seals  which  are  very  similar  to  those  which  were  empirically  studied.  However, 
any  significant  difference  in  seal  geometry  can  give  considerable  error.  Therefore, 
a  more  widely  applicable  method  has  been  sought  for  seals  of  arbitrary  geometry, 
shaft  speed,  pressure  drop,  etc. 

A  numerical  algorithm  for  solving  the  Navier-Stokes  equations  is  widely  appli¬ 
cable  for  seals.  Furthermore,  such  a  computational  tool  does  not  require  extensive 
empirical  data  as  a  user  input.  However,  this  approach  and  can  be  quite  expensive. 
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Fig.  1  Basic  configuration  and  expected  streamline 
pattern  for  a  generic  straight-through  labyrinth  seal 


StofF’s  investigation  [2]  of  incompressible  flow  in  seals  was  the  first  attempt  at  such 
a  numerical  golution.  His  work  involved  an  extension  of  the  TEACH  computer  pro- 
pram  [3],  which  constitutes  the  basic  approach  used  in  this  study  also.  Some  of  the 
more  recent  numerical  results  are  those  of  Rhode  and  Sobolik  [4]  and  Sobolik  [5] 
involving  subsonic  compressible  flow.  The  current  study  builds  upon  that  work. 

The  primary  objective  of  the  present  research  is  to  predict  the  flowfield 
development  from  cavity  to  cavity  along  a  labyrinth  seal.  Secondary  objectives 
include  comparison  of  swirl  velocity  and  other  quantities  for  corresponding  teeth- 
on-stator  and  teeth-on-rotor  generic  labyrinth  seals.  The  computer  code  of  Sobolik 
[5]  was  extended  in  various  ways  in  order  to  achieve  these  objectives. 


PREVIOUS  WORK 


Most  of  the  leakage  oriented  research  has  been  aimed  at  developing  simple 
empirical  relationships  between  the  leakage  flow  rate  and  the  overall  pressure  drop, 
se<d  geometry,  shaft  speed  and  Reynolds  number.  Details  such  as  velocity  profiles 
have  only  very  recently  become  the  objective  of  experimental  researchers. 

Work  done  on  labyrinth  seeds  can  be  classified  into  compressible  fluid  leakage, 
incompressible  fluid  leakage  and  rotordynamic  instability  studies.  The  paper  by 
Sneck  [6]  and  the  theses  by  Cogan  [7]  and  Sobilik  [5]  were  invaluable  in  compiling 
this  review.  The  topic  of  instability  is  not  addressed  in  this  paper,  and  thus  related 
studies  are  not  included. 

Compressible  Fluid  Leakage  Studies 

Martin  [8],  Stodola  [9],  Gercke  [10],  Egli  [11],  Dolin  and  Brown  [12],  Hodkinson 
[13],  Kearton  and  Keh  [14],  Zabriskie  and  Sternlicht  [15],  Vermes  [16],  Rao  and 
Naravanamurthi  [17],  Deich,  et  al.  [18]  and  Benvenuti  [19]  are  examples  of  early 
work  with  compressible  flow  seals.  These  studies  drew  conclusions  based  on  gross 
overall  characteristics  of  the  flowfield.  These  researchers  made  assumptions  such  as 
constant  flow  area,  uniform  axial  velocity  profile  or  negligible  turbulent  kinetic 
energy  carry-over.  These  assumptions  greatly  simplify  the  problem.  However, 
models  developed  in  this  way  have  been  very  limited  in  application. 

The  first  to  take  detailed  measurements  of  compressible  flow  in  labyrinth  seals 
was  Hauck  [20].  He  measured  the  axial  and  swirl  components  of  velocity  in  straight- 
through  and  stepped  seals.  Velocity  profiles  were  measured  both  in  the  cavity  as 
well  as  in  the  leakage  region.  Hauck's  test  facility  operated  at  both  4000  and  5000 
CPM,  and  allowed  for  rotor  eccentricity  up  to  75  percent  of  the  clearance  for  a 


concentric  rotor. 

Apparently  the  first  to  present  a  detailed  numerical  solution  for  the  compress¬ 
ible  flowfield  within  a  labyrinth  seal  were  Rhode  and  Sobolik  [4].  A  finite  difference 
method  was  used  employing  Patankar’s  [21]  SIMPLE  (Semi-Implicit  Method  for 
Pressure-Linked  Equations)  algorithm.  Converged  solutions  were  obtained  using 
each  of  two  different  finite  differencing  schemes:  the  Hybrid  scheme  developed  by 
Spalding  [22],  and  the  QUICK  scheme  developed  by  Leonard  [23].  Zimmermann 
and  Wolff  [24]  developed  an  analytical  model  for  correlation  of  leakage  and  pressure- 
drop.  Previously  obtained  numerical  solutions  were  used  to  integrate  cavity  exit 
velocity  profiles  and  compute  standard  loss  coefficients  for  each  cavity. 

Kirk  [25]  developed  a  relatively  simple  computer  program  to  calculate  the 
circumferential  swirl  and  pressure  distribution.  He  computed  a  nondimensionalized 
tangential  velocity  ratio  and  compared  with  previously  obtained  experimented 
results.  The  experimental  data  reveals  a  reduction  of  swirl,  as  the  flow  moves 
radially  inward,  that  was  not  predicted.  However,  predicted  and  experimental 
results  for  swirl  never  differed  by  more  than  25  percent. 

Nordmann,  et  al.  [26]  also  used  a  finite  difference  program  to  simulate  flow  in 
straight-through  labyrinth  seals.  The  effort  was  directed  toward  the  fluid  forces  and 
rotordynamic  coefficients.  But  they  also  plotted  leakage  as  a  function  of  pressure 
ratio,  comparing  against  previous  obtained  experimentally  results. 

Incompressible  Fluid  Leakage  Studies 

Jerie  [27],  Bell  and  Bergelin  [28],  Nikitin  and  Ipatov  [29]  and  Han  [30] 
are  examples  of  earlier  studies  dealing  with  gross  overall  characteristics.  The 
conclusions  reached  by  these  researchers  developed  the  application  of  flow  and 
velocity  carry-over  coefficients,  critical  Reynolds  number  for  transition  to  turbulence 


and  friction  ©efficients. 

Three  recent  studies  constitute  the  first  investigations  pertaining  to  detailed 
characteristics  of  velocity  and  pressure  throughout  a  seal  cavity.  One  was  conducted 
by  Stoff  [2].  He  used  a  finite  difference  computer  program  descended  from  the 
TEACH  program  to  solve  the  Reynolds-averaged  Navier-Stokes  equations  along 
with  those  of  the  k  —  e  turbulence  model.  Predictions  as  well  as  experimental 
measurements  of  water  in  a  large  scale  straight-through  seal  facility  were  obtained. 
He  compared  predictions  of  a  single  radial  profile  of  mean  swirl  velocity  and  rms 
swirl  velocity  with  corresponding  experimental  measurements  for  an  axial  station 
midway  between  adjacent  teeth. 

Rhode,  et  al.  [31]  is  the  second  detailed  study.  That  paper  reports  comparison 
predictions  of  two  cavity  configurations  for  a  straight-through  seal.  A  more  recent 
version  of  a  swirl-flow  finite  difference  computer  program,  also  descended  from 
TEACH,  was  employed.  Detailed  radial  profiles  of  the  three  velocity  components, 
pressure  and  turbulence  kinetic  energy  were  shown.  In  addition,  a  recent  convective 
differencing  scheme  was  evaluated  for  numerical  stability  and  accuracy. 

Finally,  Demko  [32-34]  conducted  a  study  using  both  computational  and  ex¬ 
perimental  methods.  He  presented  corresponding  predictions  and  measurements  of 
selected  quantities  such  as  axial  and  radial  velocity  and  pressure.  By  investigating 
flow  at  relatively  high  Taylor  numbers  (typically  Ta>  104)  he  was  able  to  com¬ 
putationally  and  experimentally  verify  the  existence  of  a  double  recirculation  zone 
pattern.  He  presented  his  results  in  the  form  of  a  flow  map  which  one  can  use  to 
predict  the  existence  of  either  a  single  versus  a  double  recirculation  zone  within  the 
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COMPUTATIONAL  APPROACH  , 

« 

General  Methodology 

The  compressible  flow  program  for  flow  in  labyrinth  seeds  is  a  descendant  of  the 
TEACH  program  [3].  Seven  simultaneous,  partial  differential  equations  are  solved. 
These  are:  (a)  the  compressible,  axisymmetric  form  of  the  Reynolds-averaged 
equations  for  conservation  of  momentum  (with  x,  r,  and  6  time-mean  velocity 
components  u,  v,  w),  (b)  the  two  turbulence  transport  equations  constituting  the 
k  —  e  turbulence  model,  (c)  the  energy  equation,  and  (d)  a  pressure  equation  which 
determines  the  pressure  field  and  enforces  conservation  of  mass.  Each  of  these  can 
be  divided  into  convection,  diffusion  and  source  term  components  and  cast  in  the 
general  form 


-[  — (pur^)-f  -(prrd) 


_JL{rr  JL(  r  -  <? 

dx (  T*dx )  dr^rT*dT^  5* 


(1) 


Phi  is  the  general  dependent  variable  and  is  the  diffusive  coefficient 
respective  to  each  equation.  For  example,  in  the  x-momentum  equation  <f>  —  u 
and  T  =  Me//-  The  effects  of  turbulence  are  incorporated  through  modification  of 
the  laminar-flow  momentum  diffusion  coefficient  m-  The  turbulent  viscosity  Mr  is 
added  to  m  giving  an  effective  viscosity  Me//  which  is  given  by 


Me//  =  Mr  r-  m 


(2) 


The  two-equation  k  —  e  model  f35i  was  utilized  to  evaluate  Mr- 


12 


K.  J 


>v 

V£>, 

t/v 

W 


»V-“ 

V> 

>  V 

v  V 


Figure  2  shows  the  staggered  grid  system  on  which  the  finite  difference 
equations  are  solved.  Values  for  all  variables  except  u  and  v  are  stored  at  the 
intersections  of  the  illustrated  grid  lines,  such  as  point  P.  The  axial  and  radial 
velocity  components  are  stored  at  locations  denoted  by  arrows.  Different  control 
volumes  C,  U  and  V  for  the  variables  stored  at  the  p,  w  and  s  locations  are  shown 
in  Fig.  2.  An  example  of  the  computational  domain  that  has  been  chosen  for  a 
labyrinth  seal  can  be  seen  in  Fig.  3. 

Each  4>  has  a  corresponding  finite-difference  equation  which  was  obtained  from 
Eqn.  1  by  applying  the  Gauss  Divergence  Theorem  and  expressing  the  result 
in  terms  of  neighboring  grid  point  values.  The  incorporation  of  the  differencing 
schemes  is  discussed  in  Rhode,  et  al.  [31  j.  As  an  example,  the  axial  momentum 
equation  from  which  um ,  a  velocity  estimate,  is  computed  is 


i 

To  satisfy  conservation  of  mass  locally,  the  values  for  u“  and  v'  must  undergo  a 
correction.  This  correction  is  accomplished  using  the  pressure-correction  P’  values. 
The  equation  for  this  quantity  is  derived  from  conservation  of  mass  and  momentum 
and  is  explained  in  further  detail  in  a  subsequent  section.  The  expression  for  the 
corrected  axial  velocity  formulation  is 

«,  =  «;  +  du{p'w  -  p;)  (4) 
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Du  =  — 
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The  radial  velocity  is  corrected  similarly.  The  P1  equation  takes  the  form 
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“pP'p  =  £  +  Z  (5) 

i 

as  discussed  in  Patankar  [21]. 

For  the  previously  considered  incompressible  flow  cases  [31],  the  flow  could 
be  considered  streamwise  periodic.  Thus,  cavity  inlet  boundary  conditions  for  all 
variables  excepting  pressure  were  set  equal  to  the  newly  calculated  exit  values  at 
the  conclusion  of  each  iteration.  Flow  at  the  cavity  outlet  is  strongly  convective  in 
nature,  and  thus  the  upwind  differencing  used  there  requires  no  outlet  boundary 
conditions.  For  the  present  compressible  flow  grid  independence  cases,  the  inlet 
values  resulting  from  the  above  procedure  were  adopted  in  dimensionless  form. 
Alternately,  the  cavity-by-cavity  development  computations  used  the  exit  values 
from  the  previous  successive  cavity  as  the  corresponding  inlet  values.  Further,  the 
empirically-based  law  of  the  wall  was  used  to  evaluate  normal  derivatives  of  velocity 
tangent  to  a  given  wall. 

False  Diffusion 

The  upwind  differencing  scheme  has  been  used  extensively,  however,  it  can 
introduce  false  diffusion.  Phis  truncation  error  occurs  predominantly  u:  flows  where 
convection  dominates  (i.e.,  where  the  grid  Peclet  number  |Pe!  =  exceeds  2.0) 
with  substantial  streamline-to-grid  skewness  and  substantial  diffusion  normal  to  the 
streamlines. 

False  diffusion  can  cause  an  overly  diffusive  solution.  Recirculating  flows  are 
particularly  susceptible  to  the  effects  of  false  diffusion  because  of  the  certainty 
of  considerable  velocity  gradients  and  streamline-to-grid  skewness.  The  QUICK 
(Quadratic  Upstream  Interpolation  for  Convective  Kinematics)  scheme  of  Leonard 


[23]  was  derived  with  the  intention  of  reducing  false  diffusion.  It  employs  a  three- 
point,  upwind-shifted  interpolation  formula.  Predictions  for  turbulent  flows  by 
Leschziner  and  Rodi  [36]  and  Han  et  al.  [37]  have  indicated  distinct  advantages  over 
the  Hybrid  upwind/central  differencing  scheme.  Moreover,  comparisons  between 
the  Hybrid  and  QUICK  schemes  for  incompressible  flow  in  labyrinth  seal  catties 
by  Rhode,  et  al.  [31]  for  example,  showed  that  the  QUICK  scheme  yielded  grid- 
independent  solutions  on  considerably  coarser  grids  than  for  the  Hybrid  scheme. 

QUICK  has  been  used  in  formulating  the  convective  terms  of  the  momentum 
equations  only.  Leschziner  and  Rodi  [36]  showed  the  k  and  e  solutions  are  not 
significantly  affected  by  alternative  differencing  schemes  because  of  the  source-term 
dominance  of  the  corresponding  transport  equations.  The  k  and  e  source  terms 
contain  the  large  generation  and  dissipation  effects. 

Computational  Developments 

An  improved  version  of  the  compressible-flow  P'  equation  was  developed  and 
tested.  The  essence  of  the  modifications  are  best  understood  after  reviewing  a 
derivation  of  the  P'  equation  for  a  simple  case  of  one-dimensional  axial  flow'.  The 
derivation  begins  with  the  discretized  form  of  the  steady  flow'  mass  conservation 
equation 

(Pu)e  -  (pu)u-  =  0  (6) 

Lower  case  subscripts  refer  to  values  at  control  volume  faces.  A  typical  control 
volume  for  P1  is  shown  in  Fig.  2  as  the  C  control  volume.  The  mass  flux  through 
the  east  face,  for  example,  where  the  p'eu'e  term  has  been  neglected  is 
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(^)e  =  {p*  +/>')*(«*  +»')«  =  PlK  +  PlK  +PeK  (7) 

The  quantities  p“  and  u*  are  the  latest  estimates  to  be  corrected  via  p'  and  u1  in 
order  to  satisfy  the  mass  conservation  equation. 

The  previous  version  of  the  program  used  a  formulation  for  the  p'  terms  using 
the  ideal  gas  equation  of  state  which  gives 


M  K  P'e  ^ 


(8) 


Note  that  the  expression  for  p’t  in  Eqn.  8  is  the  simple  arithmetic  mean.  It 
was  found  in  the  present  study  that  expressing  p'e  in  this  way  contributed  to 
numerical  instabilities  during  the  iterative  solution  algorithm.  This  effect  was  more 
pronounced  when  the  grid  lines  were  non-uniformly  distributed  within  the  flow 
domain.  An  alternative  formulation  of  p'e  that  reduced  numerical  instabilities  when 
using  non-uniform  grids  was  sought.  Several  variations  were  numerically  tested  for 
a  labyrinth  cavity  flowfield.  Evaluation  of  p  directly  from  the  ideal  gas  equation  of 
state  without  using  a  p'  quantity  was  an  improvement  over  the  p  =  p*  Vp'  approach. 
This  led  to  an  algorithm  that  is  less  susceptible  to  numerical  instabilities.  The 
expression  replacing  Eqn.  7  is 


(pu)e  =  Pt{u*  +  u')e  =  ptu\  +  ptu't 
here  pe  is  formulated  using  the  harmonic  mean  expression 
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(9) 


(10) 


From  the  axial  momentum  equation,  one  finds  the  following  expressions  [21] 
for  velocity  corrections 


u 


'  =  d:(p'p  -  p'E) 


(ii) 


<  =  Dl(Pw  -  P'p ) 
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Upon  appropriate  east  and  west  face  mass  flux  substitutions  into  Eqn.  6  the 
following  difference  equation  is  obtained  for  P' 


[Ae  +  Aw)Pp  =  AeP'e  +  AwP'w  ■"  Am* 


(12) 


Ae  =  PeD* 


Aw  =  pwD “ 


Am*  =  pwu'w  -  peu\ 


The  P'  equation  used  in  the  current  computer  program  was  derived  from  the 
two-dimensional  flow  mass  conservation  equation  in  cylindrical  coordinates.  The 
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procedure  exactly  parallels  that  given  here.  Appropriately,  P'  =  0  along  boundaries 
was  utilized  for  boundary  conditions. 

The  TDMA  (Tri-Diagonal  Matrix  Algorithm)  is  used  to  solve  for  the  dependent 
variable  of  each  finite  difference  equation  in  a  line-by-line  fashion.  The  previous 
version  of  the  code  solved  for  one  vertical  column  at  a  time,  beginning  at  the  west 
boundary  and  sweeping  line-by-line  to  the  east  boundary.  This  procedure  causes 
the  inlet  boundary  condition,  usually  at  the  west  boundary,  to  be  felt  quickly  by 
interior  points  of  the  calculation  domain.  This  desirable  effect  is  greatly  enhanced 
if  the  fluid  flow  is  in  the  same  direction  as  the  line-by-line  sweeping  motion  of  the 
TDMA. 

The  previous  method  of  solving  for  vertical  columns  while  sweeping  from  west 
to  east  tends  to  work  very  well  in  the  predominantly  straight-through  west  to 
east  leakage  flow  region.  It  does  not,  however,  efficiently  incorporate  the  effect 
of  boundary  conditions  within  the  recirculation  zone  of  the  cavity.  For  this  reason 
an  ADI  sweeping  procedure  was  implemented.  It  alternates  between  west-east  and 
north-south  sweeps,  solving  for  vertical  columns  and  horizontal  rows,  respectively. 
In  this  way,  information  from  the  boundary  conditions  is  efficiently  transmitted  to 
all  regions  of  the  flow  domain. 
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RESULTS  AND  DISCUSSION 


The  Seals  Considered 

In  order  to  investigate  the  cavity-by-cavity  distribution  of  various  quantities,  a 
TOS  and  also  a  TOR  seal  were  simulated.  The  TOS  cavities  had  a  shaft  radius  of 
0.0725  m,  shaft  centerline-to-cavity  base  radius  of  0.07609  m  and  shaft  centerline- 
to-tooth  periphery  radius  of  0.07291  m.  The  TOR  cavities  had  a  shaft  radius  of 
0.0725  m,  a  shaft  centerline-to-stator  radius  of  0.07609  m,  and  shaft  centerline  to 
tooth  periphery  radius  of  0.07568  m.  Both  seals  had  a  cavity  width  of  0.002825  m 
and  a  tooth  width  of  0.00035  m. 

Grid  Independence 

Figure  4  illustrates  the  effect  of  inlet  leakage  Mach  number  on  grid  indepen¬ 
dence.  A  uniform  inlet  temperature  of  293°A'  was  assumed.  After  an  inlet  leakage 
Mach  number  Mr  was  selected,  the  inlet  stator  wall  pressure  was  set  to  a  value 
yielding  an  inlet  leakage  Reynolds  number  of  Rex  =  ^2c  =  2.60  x  104.  Ax¬ 
ial  velocity  and  turbulence  kinetic  energy  profiles  were  computed  from  a  fourth- 
order  curve  fit  of  an  incompressible  flow  problem  previously  simulated  [31].  The 
shaft  rotational  speed  was  set  to  8000  CPM  and  the  radial  distribution  of  inlet 
swirl  velocity  was  taken  from  the  measurements  of  Stoff  [2],  This  gives  an  inlet 
Reg  =  2-  =  2.39  x  104  and  6.99  x  103  for  the  TOS  and  TOR  cavities,  respec¬ 

tively,  at  A/,  =  0.2.  The  respective  values  at  Mx  =  0.5  are  9.55  x  103  and  2.19  x  103. 

The  variation  of  cavity  bulk  pressure  drop  A P*  with  inlet  Mach  number  Mi  is 
shown  in  Fig.  4  for  the  teeth-on-stator  (TOS)  seal  design.  As  expected,  M{  =  0.2 
possesses  less  grid  size  dependence  than  A/,  =  0.5.  The  solution  obtained  using  a 


35x45  (x  and  r  direction)  grid  is  essentially  grid  independent  near  Mi  =  0.2.  The 
change  in  AJ at  Mi  =  0.2  was  approximately  five  percent  in  changing  from  the 
20x30  to  the  35x45  grid.  This  will  be  substantially  less  for  a  similar  comparison 
using  35x45  and  50x60  grids.  Observe  the  increased  dependency  on  grid  spacing  at 
Mi  =  0.5.  This  is  attributed  to  the  much  higher  pressure  gradients,  and  in  turn, 
density  gradients  at  the  higher  Mach  numbers. 

Almost  identical  behavior  was  found  for  the  teeth-on-rotor  (TOR)  design, 
which  is  not  shown  here.  However,  the  value  from  the  65x75  grid  at  Mi  =  0.5 
for  the  TOR  design  is  0.474,  whereas  the  corresponding  value  for  the  TOS  design 
is  0.429.  This  indicates  that  the  TOR  design  provides  around  10.5  percent  greater 
flow  resistance.  Stated  differently,  for  a  given  A P*,  the  TOR  seal  is  expected  to 
give  less  leakage.  This  agrees  with  a  corresponding  experimental  comparison  where 
both  seals  contained  fifteen  cavities  of  the  present  design. 

Figure  5  shows  the  streamline  pattern  in  a  TOS  generic  labyrinth  seal  for 
Mt  =  0.5.  It  is  similar  to  that  for  incompressible  flow’  in  a  similar  generic  seal 
in  the  paper  by  Rhode  et  al.  [31].  The  reattachment  stagnation  point  in  close 
proximity  to  the  high  speed  leakage  flow  region  gives  rise  to  large  velocity  gradients 
just  upstream  of  the  downstream  tooth. 

Cavitv-bv-Cavity  Flowfield  Development 

The  following  results  were  obtained  by  utilizing  the  finit  difference  code  to 
simulate  the  flowffeld  development  through  a  labyrinth  seal,  one  cavity  at  a  time. 
Figure  1  illustrates  the  cavity  numbering  scheme  which  is  utilized  in  the  following 
discussion.  At  the  inlet  to  the  first  cavity  a  leakage  Mach  number  of  0.3  was 
assumed.  A  uniform  temperature  profile  of  300° K  was  used  as  well  as  a  uniform 
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axial  velocity  profile.  The  inlet  radial  component,  of  velocity  was  assumed  to  be 
zero.  An  inlet  axial  leakage  Reynolds  number  of  32,000  was  assumed.  A  uniform 
inlet  profile  of  turbulence  kinetic  energy  was  used  and  a  value  of  3  percent  of  the 
mean  flow  kinetic  energy  was  specified.  The  shaft  speed  was  set  to  9500  RPM. 

Results  were  obtained  using  two  different  values  of  non-dimensional  bulk  swirl 
velocity  W’  for  comparison.  The  TOS  cases  were  run  at  W*  of  0.0  and  0.58,  while 
the  TOR  cases  had  W*  of  0.0  and  0.49.  These  non-zero  inlet  swirl  cases  give  an 
Reg  of  1.29  x  10°  and  1.14  x  105,  respectively.  These  values  of  W’  were  selected  as 
estimates  of  the  inlet  swirl  resulting  from  rotating  machinery  upstream  of  the  seal. 

After  obtaining  a  converged  solution  for  the  first  cavity,  the  field  variable  values 
at  the  cavity  exit  plane  are  used  as  inlet  values  for  the  next  cavity.  Thus  it  was 
possible  to  simulate  the  flow  development  along  the  labyrinth  seal  by  approximating 
inlet  conditions  for  the  first  cavity  and  then  proceeding  downstream  one  cavity  at 
a  time. 

For  the  TOS  case  of  zero  inlet  swirl,  Fig.  6  shows  the  increase  of  W*  with  an 
almost  constant  slope.  This  result  will  be  useful  to  the  developers  of  simple  flow 
models  for  calculating  the  net  rotordynamic  force  which  arises  when  the  rotor  is 
eccentric  with  respect  to  the  housing.  It  is  particularly  noteworthy  in  this  regard 
that  the  calculation  of  correct  rotordynamic  forces  is  sensitive  to  the  correctly 
simulated  distribution  of  W*.  Further  there  are  no  previous  measurements  or 
predictions  of  this  distribution  so  that  the  developers  of  simple  models  have  been 
unable  to  verify  this  aspect  of  their  modeling.  Thus  it  is  significant  that  Fig.  6 
shows  a  slow  rate  of  growth  of  W* .  In  fact,  such  generic  TOS  seals  with  as  many 
as  fifteen  cavities  never  even  approach  the  asymptotic  value  wThich  is  expected  to 
be  slightly  less  than  0.5.  Figure  7  shows  a  linear  decrease  for  the  TOS  case  of 
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Fig.  6  Cavity-by-cavity  development 
of  dimensionless  inlet  swirl  Wr*  for  TOS 
with  W‘  =  0.0  at  the  first  cavity 
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Fig.  7  Cavity-by-cavity  development 
of  dimensionless  inlet  swirl  W*  for  TOS 
with  IF*  =  0.58  at  the  first  cavity 
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W',*  =  0.58  and  thus  indicates  an  asymptotic  value  for  W*  which  is  slightly  less 
than  0.5. 

The  distributions  shown  in  Figs.  8  t  hough  10  reveal  changes  in  the  field 
variables  as  the  flow  proceeds  cavity-by-cavity  downstream.  The  axial  velocity 
distributions  of  Fig.  8  indicate  that  this  quantity  has  become  essentially  streamwdse 
periodic  within  the  third  cavity.  Even  though  Figs.  6  and  7  showed  an  almost  linear 
development  of  W*,  the  detailed  distributions  in  Fig.  9  show  a  significantly  non¬ 
linear  rate  of  increase  within  the  recirculation  zone  in  the  cavity  region.  The  large 
turbulent  shear  stress  Trg  on  the  shaft  gives  rise  to  the  shape  of  each  radial  profile 
near  the  shaft.  Finally,  Fig.  10  shows  the  effect  of  radial  turbulent  diffusion  on 
turbulence  energy  development  in  the  leakage  region.  This  is  most  clearly  shown 
for  the  developing  profile  in  the  first  cavity.  Although  not  shown  here,  the  radial 
velocity  profiles  showed  very  little  change. 

The  leakage  Mach  number  increased  from  a  value  of  0.30  at  the  inlet  of  the 
first  cavity  to  0.33  at  the  exit  of  the  fifth.  By  locating  Mi  =  0.33  on  Fig.  4,  one 
can  determine  that  the  35x45  grid  used  for  the  swirl  growth  investigation  was  fine 
enough  to  avoid  grid  independence  problems. 

The  TOR  results  for  the  lTt*  development  shown  in  Figs.  11  and  12  are  very 
similar  to  those  for  the  corresponding  TOS  cases.  One  important  difference  however, 
is  the  much  higher  growth  rate  for  the  TOR  case.  The  slow  linear  growth  exhibited 
in  Fig.  12  indicates  that  W?  is  close  to  its  asymptotic  value  which  from  these 
predictions  is  apparently  near  0.55.  This  result  is  somewhat  different  from  that  for 
incompressible  flow  through  a  similar  cavity  by  Rhode,  et  al.  [31]  indicating  an 
asymptotic  W’  value  of  0.65. 

Figures  13  through  15  correspond  to  Figs.  8  through  10  for  the  TOS  cases  and 
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Fig.  9  Cavity-by-cavity  dimensionless  swirl 
velocity  distribution  for  TOS  showing  the 

first  cavity ( - ),  third  cavitv( - 

and  fifth  cavity( -  -  - ) 


Fig.  10  Cavity-by-cavity  dimensionless  turbulence 
kinetic  energy  distribution  for  TOS  showing  the 

first  cavity ( - ),  third  cavity( - ) 

and  fifth  cavity( -  -  - ) 
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Fig.  12  Cavity-by-cavity  development  of 
dimensionless  inlet  swirl  W*  for  TOR 
with  IT,*  =  0.49  at  the  first  cavity 


Fig.  14  Cavity-by-cavity  dimensionless  swirl 
velocity  distribution  for  TOR  showing  the 

first  cavity( - ),  third  cavity ( - 

and  fifth  cavity( -  -  - ) 


most  of  the  results  are  nearly  identical  to  those  for  TOS.  Major  items  of  interest 
again  are  approximate  streamwise  periodicity  for  axial  velocity  at  the  third  cavity 
and  swirl  profile  development.  One  difference  between  corresponding  TOS  and  TOR 
cases  is  the  swirl  velocity  profile  in  the  leakage  flow  region.  For  the  TOS  case,  Fig. 
9  shows  a  larger  ^  than  does  Fig.  14  for  TOR.  This  is  due  to  the  fact  that  the 
large  shear  stress  rr$  on  the  shaft  acts  on  the  leakage  flow  over  the  entire  length 
from  cavity  inlet  to  outlet  for  the  TOS  case. 

Figure  16  clearly  contrasts  the  faster  swirl  growth  rate  of  TOR  over  TOS.  This 
leads  to  a  higher  asymptotic  swirl  velocity  for  TOR.  Again,  this  is  of  significance 
to  designers  interested  in  analyzing  rotordynamic  forces  on  the  shaft. 
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CLOSURE 

A  finite  difference  computer  program  for  predicting  compressible,  axisymmetric 
flow  in  labyrinth  seals  was  used  in  investigating  cavity-by-cavity  growth  of  swirl 
velocity  and  leakage  Mach  number  effect  on  grid  independence.  Due  to  a  current 
lack  of  experimental  data  from  compressible  flow  labyrinth  seals,  the  results  of  this 
computational  study  cannot  be  verified. 

The  Fanno  flow  behavior  computed  herein  for  the  compressible  flow  in  labyrinth 
seals  is  well  known.  The  steady  increase  in  Mach  number  from  one  cavity  to  the  next 
translates  into  a  significant  decrease  in  density  and  pressure.  These  in  turn  were 
found  to  give  rise  to  a  considerable  grid  dependence  sensitivity  at  Mach  numbers 
above  approximately  0.4.  Furthermore,  predictions  for  corresponding  TOR  and 
TOS  cavity  designs  wfith  the  same  flow  area  and  leakage  rate  showed  that  the  TOR 
cavity  had  approximately  ten  percent  greater  dimensionless  pressure  drop  than  the 
TOS  cavity.  Thus,  for  a  given  pressure  drop  the  TOR  cavity  is  expected  to  leak 
less. 

Several  conclusions  concerning  the  cavity-by-cavity  flowfield  development  can 
also  be  drawn  from  the  predictions  presented.  For  example,  it  was  clearly  shown 
that  the  cavity- by-cavity  swirl  velocity  development  is  fairly  slow.  It  is  anticipated 
that  these  values  will  prove  helpful  to  the  developers  of  simple  models  for  predicting 
the  rotordynamic  fluid  forces  on  such  seals  when  the  rotor  is  eccentric  relative  to 
the  stator.  Specifically,  such  previously  unavailable  values  may  serve  as  a  swirl 
development  prediction  test  case  for  such  models,  as  the  forces  are  sensitive  to 
swirl  velocity.  Moreover,  the  TOS  seals  develop  swirl  velocity  at  a  considerably 
slower  rate  than  do  TOR  seals.  This  indicates  that  TOS  seals  are  less  likely  to  have 
significant  destabilizing  fluid  forces  which  drive  the  rotor  in  the  detrimental  whirling 
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motion.  Thus,  the  greater  leakage  tendency  of  TOS  seals  offsets  this  advantage  of 
slower  swirl  velocity  development. 

In  addition,  the  prediction  of  cavity- by- cavity  flow  development  for  Mach 
numbers  near  0.3  indicates  that,  except  for  swirl  velocity,  most  quantities  are 
essentially  streamwise  periodic  after  the  second  cavity. 
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